9 #include <unordered_map>
18 #include <type_traits>
20 #ifdef EPI_DEBUG_VIRUS
28 #define EPIWORLD_VERSION_MAJOR 0
29 #define EPIWORLD_VERSION_MINOR 10
30 #define EPIWORLD_VERSION_PATCH 0
32 static const int epiworld_version_major = EPIWORLD_VERSION_MAJOR;
33 static const int epiworld_version_minor = EPIWORLD_VERSION_MINOR;
34 static const int epiworld_version_patch = EPIWORLD_VERSION_PATCH;
3820 void record_transmission(
int i,
int j,
int virus,
int i_expo_date);
3822 size_t get_n_viruses()
const;
3823 size_t get_n_tools()
const;
3825 void set_user_data(std::vector< std::string > names);
3826 void add_user_data(std::vector< epiworld_double > x);
3827 void add_user_data(epiworld_fast_uint j, epiworld_double x);
3846 MapVec_type<int,int> get_reproductive_number()
const;
3848 void get_reproductive_number(
3866 std::vector< epiworld_double > get_transition_probability(
3868 bool normalize =
true
3872 bool operator!=(
const DataBase<TSeq> & other)
const {
return !operator==(other);};
3883 void get_generation_time(
3884 std::vector< int > & agent_id,
3885 std::vector< int > & virus_id,
3886 std::vector< int > & time,
3887 std::vector< int > & gentime
3890 void get_generation_time(
5000 void record_transmission(
int i,
int j,
int virus,
int i_expo_date);
5002 size_t get_n_viruses()
const;
5003 size_t get_n_tools()
const;
5005 void set_user_data(std::vector< std::string > names);
5006 void add_user_data(std::vector< epiworld_double > x);
5007 void add_user_data(epiworld_fast_uint j, epiworld_double x);
5026 MapVec_type<int,int> get_reproductive_number()
const;
5028 void get_reproductive_number(
5046 std::vector< epiworld_double > get_transition_probability(
5048 bool normalize =
true
5052 bool operator!=(
const DataBase<TSeq> & other)
const {
return !operator==(other);};
5063 void get_generation_time(
5064 std::vector< int > & agent_id,
5065 std::vector< int > & virus_id,
5066 std::vector< int > & time,
5067 std::vector< int > & gentime
5070 void get_generation_time(
8479 epiworld_fast_uint nstates = 0u;
8481 bool verbose =
true;
8482 int current_date = 0;
8486 void dist_entities();
8488 std::chrono::time_point<std::chrono::steady_clock> time_start;
8489 std::chrono::time_point<std::chrono::steady_clock> time_end;
8492 std::chrono::duration<epiworld_double,std::micro> time_elapsed =
8493 std::chrono::duration<epiworld_double,std::micro>::zero();
8494 epiworld_fast_uint n_replicates = 0u;
8495 void chrono_start();
8498 std::vector<GlobalEvent<TSeq>> globalevents;
8501 bool use_queuing =
true;
8507 std::vector< Event<TSeq> > events = {};
8508 epiworld_fast_uint nactions = 0u;
8525 VirusPtr<TSeq> virus_,
8526 ToolPtr<TSeq> tool_,
8528 epiworld_fast_int new_state_,
8529 epiworld_fast_int queue_,
8530 EventFun<TSeq> call_,
8544 MixerFun<TSeq> susceptibility_reduction_mixer = susceptibility_reduction_mixer_default<TSeq>;
8545 MixerFun<TSeq> transmission_reduction_mixer = transmission_reduction_mixer_default<TSeq>;
8546 MixerFun<TSeq> recovery_enhancer_mixer = recovery_enhancer_mixer_default<TSeq>;
8547 MixerFun<TSeq> death_reduction_mixer = death_reduction_mixer_default<TSeq>;
8559 std::vector<epiworld_double> array_double_tmp;
8560 std::vector<Virus<TSeq> * > array_virus_tmp;
8568 virtual ~
Model() {};
8584 epiworld_double & operator()(std::string pname);
8586 size_t size()
const;
8595 void set_rand_engine(std::shared_ptr< std::mt19937 > & eng);
8596 std::shared_ptr< std::mt19937 > & get_rand_endgine();
8597 void seed(
size_t s);
8598 void set_rand_norm(epiworld_double mean, epiworld_double sd);
8599 void set_rand_unif(epiworld_double a, epiworld_double b);
8600 void set_rand_exp(epiworld_double lambda);
8601 void set_rand_gamma(epiworld_double alpha, epiworld_double beta);
8602 void set_rand_lognormal(epiworld_double mean, epiworld_double shape);
8603 void set_rand_binom(
int n, epiworld_double p);
8604 void set_rand_nbinom(
int n, epiworld_double p);
8605 void set_rand_geom(epiworld_double p);
8606 void set_rand_poiss(epiworld_double lambda);
8607 epiworld_double runif();
8608 epiworld_double runif(epiworld_double a, epiworld_double b);
8609 epiworld_double rnorm();
8610 epiworld_double rnorm(epiworld_double mean, epiworld_double sd);
8611 epiworld_double rgamma();
8612 epiworld_double rgamma(epiworld_double alpha, epiworld_double beta);
8613 epiworld_double rexp();
8614 epiworld_double rexp(epiworld_double lambda);
8615 epiworld_double rlognormal();
8616 epiworld_double rlognormal(epiworld_double mean, epiworld_double shape);
8618 int rbinom(
int n, epiworld_double p);
8620 int rnbinom(
int n, epiworld_double p);
8622 int rgeom(epiworld_double p);
8624 int rpoiss(epiworld_double lambda);
8642 void rm_virus(
size_t virus_pos);
8643 void rm_tool(
size_t tool_pos);
8644 void rm_entity(
size_t entity_id);
8663 const std::vector<int> & agents_ids,
8664 const std::vector<int> & entities_ids
8667 void load_agents_entities_ties(
8668 const int * agents_id,
8669 const int * entities_id,
8683 void agents_from_adjlist(
8687 bool directed =
false
8690 void agents_from_edgelist(
8691 const std::vector< int > & source,
8692 const std::vector< int > & target,
8697 void agents_from_adjlist(
AdjList al);
8699 bool is_directed()
const;
8711 std::vector< Entity<TSeq> > & get_entities();
8713 Entity<TSeq> & get_entity(
size_t entity_id,
int * entity_pos =
nullptr);
8716 epiworld_fast_uint n = 1000,
8717 epiworld_fast_uint k = 5,
8719 epiworld_double p = .01
8721 void agents_empty_graph(epiworld_fast_uint n = 1000);
8734 void update_state();
8735 void mutate_virus();
8738 epiworld_fast_uint ndays,
8742 epiworld_fast_uint ndays,
8743 epiworld_fast_uint nexperiments,
8745 std::function<
void(
size_t,
Model<TSeq>*)> fun = make_save_run<TSeq>(),
8747 bool verbose =
true,
8754 epiworld_fast_uint get_ndays()
const;
8755 epiworld_fast_uint get_n_replicates()
const;
8756 size_t get_n_entities()
const;
8757 void set_ndays(epiworld_fast_uint ndays);
8758 bool get_verbose()
const;
8776 void set_rewire_prop(epiworld_double prop);
8777 epiworld_double get_rewire_prop()
const;
8794 std::string fn_virus_info,
8795 std::string fn_virus_hist,
8796 std::string fn_tool_info,
8797 std::string fn_tool_hist,
8798 std::string fn_total_hist,
8799 std::string fn_transmission,
8800 std::string fn_transition,
8801 std::string fn_reproductive_number,
8802 std::string fn_generation_time
8816 void write_edgelist(
8820 void write_edgelist(
8821 std::vector< int > & source,
8822 std::vector< int > & target
8826 std::map<std::string, epiworld_double> & params();
8840 const Model<TSeq> & print(
bool lite =
false)
const;
8857 void add_state(std::string lab, UpdateFun<TSeq> fun =
nullptr);
8858 const std::vector< std::string > & get_states()
const;
8859 size_t get_n_states()
const;
8860 const std::vector< UpdateFun<TSeq> > & get_state_fun()
const;
8861 void print_state_codes()
const;
8873 std::vector< double > ,
8912 epiworld_double add_param(
8913 epiworld_double initial_val, std::string pname,
bool overwrite =
false
8915 Model<TSeq> & read_params(std::string fn,
bool overwrite =
false);
8916 epiworld_double get_param(epiworld_fast_uint k);
8917 epiworld_double get_param(std::string pname);
8919 void set_param(std::string pname, epiworld_double val);
8921 epiworld_double par(std::string pname)
const;
8925 std::string unit =
"auto",
8926 epiworld_double * last_elapsed =
nullptr,
8927 epiworld_double * total_elapsed =
nullptr,
8928 std::string * unit_abbr =
nullptr,
8939 void add_user_data(epiworld_fast_uint j, epiworld_double x);
8940 void add_user_data(std::vector< epiworld_double > x);
8957 std::string name =
"A global action",
8961 void add_globalevent(
8971 void run_globalevents();
8973 void clear_state_set();
8996 void set_susceptibility_reduction_mixer(MixerFun<TSeq> fun);
8997 void set_transmission_reduction_mixer(MixerFun<TSeq> fun);
8998 void set_recovery_enhancer_mixer(MixerFun<TSeq> fun);
8999 void set_death_reduction_mixer(MixerFun<TSeq> fun);
9002 const std::vector< VirusPtr<TSeq> > & get_viruses()
const;
9003 const std::vector< ToolPtr<TSeq> > & get_tools()
const;
9019 double * get_agents_data();
9020 size_t get_agents_data_ncols()
const;
9028 std::string get_name()
const;
9031 bool operator!=(
const Model<TSeq> & other)
const {
return !operator==(other);};
9048 DiagramType diagram_type = DiagramType::Mermaid,
9049 const std::string & fn_output =
"",
9795 epiworld_fast_uint nstates = 0u;
9797 bool verbose =
true;
9798 int current_date = 0;
9802 void dist_entities();
9804 std::chrono::time_point<std::chrono::steady_clock> time_start;
9805 std::chrono::time_point<std::chrono::steady_clock> time_end;
9808 std::chrono::duration<epiworld_double,std::micro> time_elapsed =
9809 std::chrono::duration<epiworld_double,std::micro>::zero();
9810 epiworld_fast_uint n_replicates = 0u;
9811 void chrono_start();
9814 std::vector<GlobalEvent<TSeq>> globalevents;
9817 bool use_queuing =
true;
9823 std::vector< Event<TSeq> > events = {};
9824 epiworld_fast_uint nactions = 0u;
9841 VirusPtr<TSeq> virus_,
9842 ToolPtr<TSeq> tool_,
9844 epiworld_fast_int new_state_,
9845 epiworld_fast_int queue_,
9846 EventFun<TSeq> call_,
9860 MixerFun<TSeq> susceptibility_reduction_mixer = susceptibility_reduction_mixer_default<TSeq>;
9861 MixerFun<TSeq> transmission_reduction_mixer = transmission_reduction_mixer_default<TSeq>;
9862 MixerFun<TSeq> recovery_enhancer_mixer = recovery_enhancer_mixer_default<TSeq>;
9863 MixerFun<TSeq> death_reduction_mixer = death_reduction_mixer_default<TSeq>;
9875 std::vector<epiworld_double> array_double_tmp;
9876 std::vector<Virus<TSeq> * > array_virus_tmp;
9884 virtual ~
Model() {};
9900 epiworld_double & operator()(std::string pname);
9902 size_t size()
const;
9911 void set_rand_engine(std::shared_ptr< std::mt19937 > & eng);
9912 std::shared_ptr< std::mt19937 > & get_rand_endgine();
9913 void seed(
size_t s);
9914 void set_rand_norm(epiworld_double mean, epiworld_double sd);
9915 void set_rand_unif(epiworld_double a, epiworld_double b);
9916 void set_rand_exp(epiworld_double lambda);
9917 void set_rand_gamma(epiworld_double alpha, epiworld_double beta);
9918 void set_rand_lognormal(epiworld_double mean, epiworld_double shape);
9919 void set_rand_binom(
int n, epiworld_double p);
9920 void set_rand_nbinom(
int n, epiworld_double p);
9921 void set_rand_geom(epiworld_double p);
9922 void set_rand_poiss(epiworld_double lambda);
9923 epiworld_double runif();
9924 epiworld_double runif(epiworld_double a, epiworld_double b);
9925 epiworld_double rnorm();
9926 epiworld_double rnorm(epiworld_double mean, epiworld_double sd);
9927 epiworld_double rgamma();
9928 epiworld_double rgamma(epiworld_double alpha, epiworld_double beta);
9929 epiworld_double rexp();
9930 epiworld_double rexp(epiworld_double lambda);
9931 epiworld_double rlognormal();
9932 epiworld_double rlognormal(epiworld_double mean, epiworld_double shape);
9934 int rbinom(
int n, epiworld_double p);
9936 int rnbinom(
int n, epiworld_double p);
9938 int rgeom(epiworld_double p);
9940 int rpoiss(epiworld_double lambda);
9958 void rm_virus(
size_t virus_pos);
9959 void rm_tool(
size_t tool_pos);
9960 void rm_entity(
size_t entity_id);
9973 void load_agents_entities_ties(std::string fn,
int skip);
9978 void load_agents_entities_ties(
9979 const std::vector<int> & agents_ids,
9980 const std::vector<int> & entities_ids
9983 void load_agents_entities_ties(
9984 const int * agents_id,
9985 const int * entities_id,
9999 void agents_from_adjlist(
10003 bool directed =
false
10006 void agents_from_edgelist(
10007 const std::vector< int > & source,
10008 const std::vector< int > & target,
10013 void agents_from_adjlist(
AdjList al);
10015 bool is_directed()
const;
10017 std::vector< Agent<TSeq> > & get_agents();
10021 std::vector< epiworld_fast_uint > get_agents_states()
const;
10023 std::vector< Viruses_const<TSeq> > get_agents_viruses()
const;
10025 std::vector< Viruses<TSeq> > get_agents_viruses();
10027 std::vector< Entity<TSeq> > & get_entities();
10029 Entity<TSeq> & get_entity(
size_t entity_id,
int * entity_pos =
nullptr);
10032 epiworld_fast_uint n = 1000,
10033 epiworld_fast_uint k = 5,
10035 epiworld_double p = .01
10037 void agents_empty_graph(epiworld_fast_uint n = 1000);
10050 void update_state();
10051 void mutate_virus();
10054 epiworld_fast_uint ndays,
10058 epiworld_fast_uint ndays,
10059 epiworld_fast_uint nexperiments,
10061 std::function<
void(
size_t,
Model<TSeq>*)> fun = make_save_run<TSeq>(),
10063 bool verbose =
true,
10068 size_t get_n_viruses()
const;
10069 size_t get_n_tools()
const;
10070 epiworld_fast_uint get_ndays()
const;
10071 epiworld_fast_uint get_n_replicates()
const;
10072 size_t get_n_entities()
const;
10073 void set_ndays(epiworld_fast_uint ndays);
10074 bool get_verbose()
const;
10092 void set_rewire_prop(epiworld_double prop);
10093 epiworld_double get_rewire_prop()
const;
10110 std::string fn_virus_info,
10111 std::string fn_virus_hist,
10112 std::string fn_tool_info,
10113 std::string fn_tool_hist,
10114 std::string fn_total_hist,
10115 std::string fn_transmission,
10116 std::string fn_transition,
10117 std::string fn_reproductive_number,
10118 std::string fn_generation_time
10132 void write_edgelist(
10136 void write_edgelist(
10137 std::vector< int > & source,
10138 std::vector< int > & target
10142 std::map<std::string, epiworld_double> & params();
10155 virtual void reset();
10156 const Model<TSeq> & print(
bool lite =
false)
const;
10173 void add_state(std::string lab, UpdateFun<TSeq> fun =
nullptr);
10174 const std::vector< std::string > & get_states()
const;
10175 size_t get_n_states()
const;
10176 const std::vector< UpdateFun<TSeq> > & get_state_fun()
const;
10177 void print_state_codes()
const;
10189 std::vector< double > ,
10228 epiworld_double add_param(
10229 epiworld_double initial_val, std::string pname,
bool overwrite =
false
10231 Model<TSeq> & read_params(std::string fn,
bool overwrite =
false);
10232 epiworld_double get_param(epiworld_fast_uint k);
10233 epiworld_double get_param(std::string pname);
10235 void set_param(std::string pname, epiworld_double val);
10237 epiworld_double par(std::string pname)
const;
10241 std::string unit =
"auto",
10242 epiworld_double * last_elapsed =
nullptr,
10243 epiworld_double * total_elapsed =
nullptr,
10244 std::string * unit_abbr =
nullptr,
10254 void set_user_data(std::vector< std::string > names);
10255 void add_user_data(epiworld_fast_uint j, epiworld_double x);
10256 void add_user_data(std::vector< epiworld_double > x);
10271 void add_globalevent(
10273 std::string name =
"A global action",
10277 void add_globalevent(
10284 void rm_globalevent(std::string name);
10285 void rm_globalevent(
size_t i);
10287 void run_globalevents();
10289 void clear_state_set();
10301 bool is_queuing_on()
const;
10312 void set_susceptibility_reduction_mixer(MixerFun<TSeq> fun);
10313 void set_transmission_reduction_mixer(MixerFun<TSeq> fun);
10314 void set_recovery_enhancer_mixer(MixerFun<TSeq> fun);
10315 void set_death_reduction_mixer(MixerFun<TSeq> fun);
10318 const std::vector< VirusPtr<TSeq> > & get_viruses()
const;
10319 const std::vector< ToolPtr<TSeq> > & get_tools()
const;
10334 void set_agents_data(
double * data_,
size_t ncols_);
10335 double * get_agents_data();
10336 size_t get_agents_data_ncols()
const;
10343 void set_name(std::string name);
10344 std::string get_name()
const;
10346 bool operator==(
const Model<TSeq> & other)
const;
10347 bool operator!=(
const Model<TSeq> & other)
const {
return !operator==(other);};
10364 DiagramType diagram_type = DiagramType::Mermaid,
10365 const std::string & fn_output =
"",
12219 template<
typename TSeq>
12227 std::ifstream filei(fn);
12230 throw std::logic_error(
"The file " + fn +
" was not found.");
12233 std::vector< std::vector< epiworld_fast_uint > > target_(entities.size());
12235 target_.reserve(1e5);
12237 while (!filei.eof())
12240 if (linenum++ < skip)
12247 throw std::logic_error(
12248 "I/O error while reading the file " +
12255 if (i >=
static_cast<int>(this->size()))
12256 throw std::range_error(
12257 "The agent["+std::to_string(linenum)+
"] = " + std::to_string(i) +
12258 " is above the max id " + std::to_string(this->size() - 1)
12261 if (j >=
static_cast<int>(this->entities.size()))
12262 throw std::range_error(
12263 "The entity["+std::to_string(linenum)+
"] = " + std::to_string(j) +
12264 " is above the max id " + std::to_string(this->entities.size() - 1)
12267 target_[j].push_back(i);
12269 population[i].add_entity(entities[j],
nullptr);
12277 template<
typename TSeq>
12279 const std::vector< int > & agents_ids,
12280 const std::vector< int > & entities_ids
12284 if (agents_ids.size() != entities_ids.size())
12285 throw std::length_error(
12286 std::string(
"The size of agents_ids (") +
12287 std::to_string(agents_ids.size()) +
12288 std::string(
") and entities_ids (") +
12289 std::to_string(entities_ids.size()) +
12290 std::string(
") must be the same.")
12293 return this->load_agents_entities_ties(
12295 entities_ids.data(),
12301 template<
typename TSeq>
12303 const int * agents_ids,
12304 const int * entities_ids,
12308 auto get_agent = [agents_ids](
int i) ->
int {
12309 return *(agents_ids + i);
12312 auto get_entity = [entities_ids](
int i) ->
int {
12313 return *(entities_ids + i);
12316 for (
size_t i = 0u; i < n; ++i)
12319 if (get_agent(i) < 0)
12320 throw std::length_error(
12321 std::string(
"agents_ids[") +
12322 std::to_string(i) +
12323 std::string(
"] = ") +
12324 std::to_string(get_agent(i)) +
12325 std::string(
" is negative.")
12328 if (get_entity(i) < 0)
12329 throw std::length_error(
12330 std::string(
"entities_ids[") +
12331 std::to_string(i) +
12332 std::string(
"] = ") +
12333 std::to_string(get_entity(i)) +
12334 std::string(
" is negative.")
12337 int pop_size =
static_cast<int>(this->population.size());
12338 if (get_agent(i) >= pop_size)
12339 throw std::length_error(
12340 std::string(
"agents_ids[") +
12341 std::to_string(i) +
12342 std::string(
"] = ") +
12343 std::to_string(get_agent(i)) +
12344 std::string(
" is out of range (population size: ") +
12345 std::to_string(pop_size) +
12349 int ent_size =
static_cast<int>(this->entities.size());
12350 if (get_entity(i) >= ent_size)
12351 throw std::length_error(
12352 std::string(
"entities_ids[") +
12353 std::to_string(i) +
12354 std::string(
"] = ") +
12355 std::to_string(get_entity(i)) +
12356 std::string(
" is out of range (entities size: ") +
12357 std::to_string(ent_size) +
12362 this->population[get_agent(i)].add_entity(
12363 this->entities[get_entity(i)],
12374 template<
typename TSeq>
12384 this->agents_from_adjlist(al);
12388 template<
typename TSeq>
12390 const std::vector< int > & source,
12391 const std::vector< int > & target,
12396 AdjList al(source, target, size, directed);
12397 agents_from_adjlist(al);
12401 template<
typename TSeq>
12405 agents_empty_graph(al.
vcount());
12407 const auto & tmpdat = al.get_dat();
12409 for (
size_t i = 0u; i < tmpdat.size(); ++i)
12413 population[i].model =
this;
12415 for (
const auto & link: tmpdat[i])
12418 population[i].add_neighbor(
12419 population[link.first],
12428 for (
auto & p: population)
12430 if (p.id >=
static_cast<int>(al.
vcount()))
12431 throw std::logic_error(
12432 "Agent's id cannot be negative above or equal to the number of agents!");
12438 template<
typename TSeq>
12441 if (population.size() == 0u)
12442 throw std::logic_error(
"The population hasn't been initialized.");
12447 template<
typename TSeq>
12453 return this->current_date;
12456 template<
typename TSeq>
12461 for (
auto & p : population)
12463 if ((p.state >= nstates) || (p.state < 0))
12465 throw std::range_error(
12466 "The agent " + std::to_string(p.id) +
12467 " has state " + std::to_string(p.state) +
12468 " which is above the maximum state of " +
12469 std::to_string(nstates - 1) +
"."
12473 if ((p.state_prev >= nstates) || (p.state_prev < 0))
12475 throw std::range_error(
12476 "The agent " + std::to_string(p.id) +
12477 " has previous state " + std::to_string(p.state_prev) +
12478 " which is above the maximum state of " +
12479 std::to_string(nstates - 1) +
"."
12488 ++this->current_date;
12491 if ((this->current_date >= 1) && verbose)
12497 template<
typename TSeq>
12499 epiworld_fast_uint ndays,
12505 throw std::logic_error(
"There are no agents in this model!");
12508 throw std::logic_error(
12509 std::string(
"No states registered in this model. ") +
12510 std::string(
"At least one state should be included. See the function -Model::add_state()-")
12514 this->ndays = ndays;
12517 engine->seed(seed);
12519 array_double_tmp.resize(std::max(
12521 static_cast<size_t>(1024)
12525 array_virus_tmp.resize(1024);
12529 epiworld_fast_int _init, _end, _removed;
12530 int nstate_int =
static_cast<int>(nstates);
12531 for (
auto & v : viruses)
12533 v->get_state(&_init, &_end, &_removed);
12536 if (((_init != -99) && (_init < 0)) || (_init >= nstate_int))
12537 throw std::range_error(
"States must be between 0 and " +
12538 std::to_string(nstates - 1));
12541 if (((_end != -99) && (_end < 0)) || (_end >= nstate_int))
12542 throw std::range_error(
"States must be between 0 and " +
12543 std::to_string(nstates - 1));
12545 if (((_removed != -99) && (_removed < 0)) || (_removed >= nstate_int))
12546 throw std::range_error(
"States must be between 0 and " +
12547 std::to_string(nstates - 1));
12551 for (
auto & t : tools)
12553 t->get_state(&_init, &_end);
12556 if (((_init != -99) && (_init < 0)) || (_init >= nstate_int))
12557 throw std::range_error(
"States must be between 0 and " +
12558 std::to_string(nstates - 1));
12561 if (((_end != -99) && (_end < 0)) || (_end >= nstate_int))
12562 throw std::range_error(
"States must be between 0 and " +
12563 std::to_string(nstates - 1));
12572 EPIWORLD_RUN((*
this))
12576 db.n_transmissions_potential = 0;
12577 db.n_transmissions_today = 0;
12582 this->update_state();
12585 this->run_globalevents();
12595 this->mutate_virus();
12600 this->current_date--;
12608 template<
typename TSeq>
12610 epiworld_fast_uint ndays,
12611 epiworld_fast_uint nexperiments,
12628 std::vector< int > seeds_n(nexperiments);
12629 for (
auto & s : seeds_n)
12631 s =
static_cast<int>(
12633 runif() *
static_cast<double>(std::numeric_limits<int>::max())
12641 EPI_DEBUG_NOTIFY_ACTIVE()
12644 bool old_verb = this->verbose;
12653 omp_set_num_threads(nthreads);
12657 static_cast<size_t>(nthreads) > nexperiments ? nexperiments : nthreads;
12660 std::vector< Model<TSeq> * > these(
12661 std::max(nthreads - 1, 0)
12664 #pragma omp parallel for shared(these, nthreads)
12665 for (
size_t i = 0u; i < static_cast<size_t>(nthreads); ++i)
12671 these[i - 1] = clone_ptr();
12677 std::vector< size_t > nreplicates(nthreads, 0);
12678 std::vector< size_t > nreplicates_csum(nthreads, 0);
12680 size_t base_replicates = nexperiments / nthreads;
12681 size_t remainder = nexperiments % nthreads;
12684 for (
int i = 0; i < nthreads; ++i)
12687 nreplicates[i] = base_replicates + (
static_cast<size_t>(i) < remainder ? 1 : 0);
12690 nreplicates_csum[i] = sums;
12691 sums += nreplicates[i];
12696 EPIWORLD_PROGRESS_BAR_WIDTH
12703 "Starting multiple runs (%i) using %i thread(s)\n",
12704 static_cast<int>(nexperiments),
12705 static_cast<int>(nthreads)
12708 pb_multiple.start();
12715 for (
size_t i = 1; i < static_cast<size_t>(std::max(nthreads - 1, 0)); ++i)
12718 if (db != these[i]->db)
12720 throw std::runtime_error(
12721 "The initial state of the models is not the same"
12727 #pragma omp parallel shared(these, nreplicates, nreplicates_csum, seeds_n) \
12728 firstprivate(nexperiments, nthreads, fun, reset, verbose, pb_multiple, ndays) \
12732 auto iam = omp_get_thread_num();
12734 for (
size_t n = 0u; n < nreplicates[iam]; ++n)
12736 size_t sim_id = nreplicates_csum[iam] + n;
12741 EPI_CHECK_USER_INTERRUPT(n);
12744 run(ndays, seeds_n[sim_id]);
12751 pb_multiple.next();
12756 these[iam - 1]->run(ndays, seeds_n[sim_id]);
12759 fun(sim_id, these[iam - 1]);
12768 n_replicates += (nexperiments - nreplicates[0u]);
12770 #pragma omp parallel for shared(these)
12771 for (
int i = 1; i < nthreads; ++i)
12773 delete these[i - 1];
12782 EPIWORLD_PROGRESS_BAR_WIDTH
12789 "Starting multiple runs (%i)\n",
12790 static_cast<int>(nexperiments)
12793 pb_multiple.start();
12797 for (
size_t n = 0u; n < nexperiments; ++n)
12801 EPI_CHECK_USER_INTERRUPT(n);
12803 run(ndays, seeds_n[n]);
12809 pb_multiple.next();
12821 template<
typename TSeq>
12828 for (
auto & p: population)
12829 if (queue[++i] > 0)
12831 if (state_fun[p.state])
12832 state_fun[p.state](&p,
this);
12839 for (
auto & p: population)
12840 if (state_fun[p.state])
12841 state_fun[p.state](&p,
this);
12849 template<
typename TSeq>
12853 size_t nmutates = 0u;
12854 for (
const auto & v: viruses)
12855 if (v->virus_functions->mutation)
12858 if (nmutates == 0u)
12865 for (
auto & p: population)
12868 if (queue[++i] == 0)
12871 if (p.virus !=
nullptr)
12872 p.virus->mutate(
this);
12880 for (
auto & p: population)
12883 if (p.virus !=
nullptr)
12884 p.virus->mutate(
this);
12893 template<
typename TSeq>
12898 template<
typename TSeq>
12900 return tools.size();
12903 template<
typename TSeq>
12908 template<
typename TSeq>
12911 return n_replicates;
12914 template<
typename TSeq>
12916 return entities.size();
12919 template<
typename TSeq>
12921 this->ndays = ndays;
12924 template<
typename TSeq>
12929 template<
typename TSeq>
12935 template<
typename TSeq>
12941 template<
typename TSeq>
12948 template<
typename TSeq>
12953 throw std::range_error(
"Proportions cannot be negative.");
12956 throw std::range_error(
"Proportions cannot be above 1.0.");
12958 rewire_prop = prop;
12961 template<
typename TSeq>
12963 return rewire_prop;
12966 template<
typename TSeq>
12970 rewire_fun(&population,
this, rewire_prop);
12974 template<
typename TSeq>
12976 std::string fn_virus_info,
12977 std::string fn_virus_hist,
12978 std::string fn_tool_info,
12979 std::string fn_tool_hist,
12980 std::string fn_total_hist,
12981 std::string fn_transmission,
12982 std::string fn_transition,
12983 std::string fn_reproductive_number,
12984 std::string fn_generation_time
12989 fn_virus_info, fn_virus_hist,
12990 fn_tool_info, fn_tool_hist,
12991 fn_total_hist, fn_transmission, fn_transition,
12992 fn_reproductive_number, fn_generation_time
12997 template<
typename TSeq>
13004 std::vector< const Agent<TSeq> * > wseq(size());
13005 for (
const auto & p: population)
13008 std::ofstream efile(fn, std::ios_base::out);
13009 efile <<
"source target\n";
13010 if (this->is_directed())
13013 for (
const auto & p : wseq)
13016 if (p->neighbors ==
nullptr)
13019 for (
auto & n : *p->neighbors)
13020 efile << p->id <<
" " << n <<
"\n";
13025 for (
const auto & p : wseq)
13028 if (p->neighbors ==
nullptr)
13031 for (
auto & n : *p->neighbors)
13032 if (
static_cast<int>(p->id) <=
static_cast<int>(n))
13033 efile << p->id <<
" " << n <<
"\n";
13040 template<
typename TSeq>
13042 std::vector< int > & source,
13043 std::vector< int > & target
13047 std::vector< const Agent<TSeq> * > wseq(size());
13048 for (
const auto & p: population)
13051 if (this->is_directed())
13054 for (
const auto & p : wseq)
13056 if (p->neighbors ==
nullptr)
13059 for (
auto & n : *p->neighbors)
13061 source.push_back(
static_cast<int>(p->id));
13062 target.push_back(
static_cast<int>(n));
13068 for (
const auto & p : wseq)
13071 if (p->neighbors ==
nullptr)
13074 for (
auto & n : *p->neighbors) {
13075 if (
static_cast<int>(p->id) <=
static_cast<int>(n)) {
13076 source.push_back(
static_cast<int>(p->id));
13077 target.push_back(
static_cast<int>(n));
13087 template<
typename TSeq>
13093 template<
typename TSeq>
13099 if (population_backup.size())
13101 population = population_backup;
13104 for (
auto & p : population)
13108 for (
size_t i = 0; i < population.size(); ++i)
13111 if (population[i] != (population_backup)[i])
13112 throw std::logic_error(
"Model::reset population doesn't match.");
13119 for (
auto & p : population)
13123 for (
auto & a: population)
13125 if (a.get_state() != 0u)
13126 throw std::logic_error(
"Model::reset population doesn't match."
13127 "Some agents are not in the baseline state.");
13131 if (entities_backup.size())
13133 entities = entities_backup;
13136 for (
size_t i = 0; i < entities.size(); ++i)
13139 if (entities[i] != (entities_backup)[i])
13140 throw std::logic_error(
"Model::reset entities don't match.");
13147 for (
auto & e: entities)
13164 initial_states_fun(
this);
17977 m->array_double_tmp[nviruses_tmp] =
17978 (1.0 - p->get_susceptibility_reduction(v, m)) *
17979 v->get_prob_infecting(m) *
17980 (1.0 - neighbor->get_transmission_reduction(v, m))
17983 m->array_virus_tmp[nviruses_tmp++] = &(*v);
17988 if (nviruses_tmp == 0u)
17992 int which = roulette(nviruses_tmp, m);
17997 p->set_virus(*m->array_virus_tmp[which], m);
18007 std::shared_ptr<std::vector<bool>> exclude_agent_bool =
18008 std::make_shared<std::vector<bool>>(0);
18010 std::shared_ptr<std::vector<epiworld_fast_uint>> exclude_agent_bool_idx =
18011 std::make_shared<std::vector<epiworld_fast_uint>>(exclude);
18018 if (exclude_agent_bool->size() == 0u)
18021 exclude_agent_bool->resize(m->get_states().size(),
false);
18022 for (
auto s : *exclude_agent_bool_idx)
18024 if (s >= exclude_agent_bool->size())
18025 throw std::logic_error(
18026 std::string(
"You are trying to exclude a state that is out of range: ") +
18027 std::to_string(s) + std::string(
". There are only ") +
18028 std::to_string(exclude_agent_bool->size()) +
18029 std::string(
" states in the model.")
18032 exclude_agent_bool->operator[](s) =
true;
18038 if (p->get_virus() !=
nullptr)
18039 throw std::logic_error(
18040 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense! ") +
18041 std::string(
"Agent id ") + std::to_string(p->get_id()) +
18042 std::string(
" has a virus.")
18046 size_t nviruses_tmp = 0u;
18047 for (
auto & neighbor: p->get_neighbors())
18051 if (exclude_agent_bool->operator[](neighbor->get_state()))
18054 auto & v = neighbor->get_virus();
18060 m->array_double_tmp[nviruses_tmp] =
18061 (1.0 - p->get_susceptibility_reduction(v, m)) *
18062 v->get_prob_infecting(m) *
18063 (1.0 - neighbor->get_transmission_reduction(v, m))
18066 m->array_virus_tmp[nviruses_tmp++] = &(*v);
18071 if (nviruses_tmp == 0u)
18075 int which = roulette(nviruses_tmp, m);
18080 p->set_virus(*m->array_virus_tmp[which], m);
18105 template<
typename TSeq = EPI_DEFAULT_TSEQ>
18107 std::vector< epiworld_fast_uint > exclude = {}
18110 if (exclude.size() == 0u)
18116 if (p->get_virus() !=
nullptr)
18117 throw std::logic_error(
18118 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense! ") +
18119 std::string(
"Agent id ") + std::to_string(p->get_id()) +
18120 std::string(
" has a virus.")
18124 size_t nviruses_tmp = 0u;
18125 for (
auto & neighbor: p->get_neighbors())
18128 if (neighbor->get_virus() ==
nullptr)
18131 auto & v = neighbor->get_virus();
18134 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
18135 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
18139 m->array_double_tmp[nviruses_tmp] =
18140 (1.0 - p->get_susceptibility_reduction(v, m)) *
18141 v->get_prob_infecting(m) *
18142 (1.0 - neighbor->get_transmission_reduction(v, m))
18145 m->array_virus_tmp[nviruses_tmp++] = &(*v);
18150 if (nviruses_tmp == 0u)
18154 int which = roulette(nviruses_tmp, m);
18159 return m->array_virus_tmp[which];
18169 std::shared_ptr<std::vector<bool>> exclude_agent_bool =
18170 std::make_shared<std::vector<bool>>(0);
18172 std::shared_ptr<std::vector<epiworld_fast_uint>> exclude_agent_bool_idx =
18173 std::make_shared<std::vector<epiworld_fast_uint>>(exclude);
18180 if (exclude_agent_bool->size() == 0u)
18183 exclude_agent_bool->resize(m->get_states().size(),
false);
18184 for (
auto s : *exclude_agent_bool_idx)
18186 if (s >= exclude_agent_bool->size())
18187 throw std::logic_error(
18188 std::string(
"You are trying to exclude a state that is out of range: ") +
18189 std::to_string(s) + std::string(
". There are only ") +
18190 std::to_string(exclude_agent_bool->size()) +
18191 std::string(
" states in the model.")
18194 exclude_agent_bool->operator[](s) =
true;
18200 if (p->get_virus() !=
nullptr)
18201 throw std::logic_error(
18202 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense! ") +
18203 std::string(
"Agent id ") + std::to_string(p->get_id()) +
18204 std::string(
" has a virus.")
18208 size_t nviruses_tmp = 0u;
18209 for (
auto & neighbor: p->get_neighbors())
18213 if (exclude_agent_bool->operator[](neighbor->get_state()))
18216 if (neighbor->get_virus() ==
nullptr)
18219 auto & v = neighbor->get_virus();
18222 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
18223 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
18227 m->array_double_tmp[nviruses_tmp] =
18228 (1.0 - p->get_susceptibility_reduction(v, m)) *
18229 v->get_prob_infecting(m) *
18230 (1.0 - neighbor->get_transmission_reduction(v, m))
18233 m->array_virus_tmp[nviruses_tmp++] = &(*v);
18238 if (nviruses_tmp == 0u)
18242 int which = roulette(nviruses_tmp, m);
18247 return m->array_virus_tmp[which];
18273 template<
typename TSeq = EPI_DEFAULT_TSEQ>
18277 if (p->get_virus() !=
nullptr)
18278 throw std::logic_error(
18279 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense!") +
18280 std::string(
"Agent id ") + std::to_string(p->
get_id()) +
18281 std::string(
" has a virus.")
18285 size_t nviruses_tmp = 0u;
18286 for (
auto & neighbor: p->get_neighbors())
18289 int _vcount_neigh = 0;
18292 if (neighbor->get_virus() ==
nullptr)
18295 auto & v = neighbor->get_virus();
18298 if (nviruses_tmp >= m->array_virus_tmp.size())
18299 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
18303 m->array_double_tmp[nviruses_tmp] =
18304 (1.0 - p->get_susceptibility_reduction(v, m)) *
18305 v->get_prob_infecting(m) *
18306 (1.0 - neighbor->get_transmission_reduction(v, m))
18309 m->array_virus_tmp[nviruses_tmp++] = &(*v);
18313 (m->array_double_tmp[nviruses_tmp - 1] < 0.0) |
18314 (m->array_double_tmp[nviruses_tmp - 1] > 1.0)
18318 "[epi-debug] Agent %i's virus %i has transmission prob outside of [0, 1]: %.4f!\n",
18319 static_cast<int>(neighbor->get_id()),
18320 static_cast<int>(_vcount_neigh++),
18321 m->array_double_tmp[nviruses_tmp - 1]
18330 if (nviruses_tmp == 0u)
18334 m->get_db().n_transmissions_potential++;
18338 int which = roulette(nviruses_tmp, m);
18344 m->get_db().n_transmissions_today++;
18347 return m->array_virus_tmp[which];
18561 epiworld_fast_uint nstates = 0u;
18563 bool verbose =
true;
18564 int current_date = 0;
18568 void dist_entities();
18570 std::chrono::time_point<std::chrono::steady_clock> time_start;
18571 std::chrono::time_point<std::chrono::steady_clock> time_end;
18574 std::chrono::duration<epiworld_double,std::micro> time_elapsed =
18575 std::chrono::duration<epiworld_double,std::micro>::zero();
18576 epiworld_fast_uint n_replicates = 0u;
18577 void chrono_start();
18580 std::vector<GlobalEvent<TSeq>> globalevents;
18583 bool use_queuing =
true;
18589 std::vector< Event<TSeq> > events = {};
18590 epiworld_fast_uint nactions = 0u;
18607 VirusPtr<TSeq> virus_,
18608 ToolPtr<TSeq> tool_,
18610 epiworld_fast_int new_state_,
18611 epiworld_fast_int queue_,
18612 EventFun<TSeq> call_,
18626 MixerFun<TSeq> susceptibility_reduction_mixer = susceptibility_reduction_mixer_default<TSeq>;
18627 MixerFun<TSeq> transmission_reduction_mixer = transmission_reduction_mixer_default<TSeq>;
18628 MixerFun<TSeq> recovery_enhancer_mixer = recovery_enhancer_mixer_default<TSeq>;
18629 MixerFun<TSeq> death_reduction_mixer = death_reduction_mixer_default<TSeq>;
18641 std::vector<epiworld_double> array_double_tmp;
18642 std::vector<Virus<TSeq> * > array_virus_tmp;
18650 virtual ~
Model() {};
18666 epiworld_double & operator()(std::string pname);
18668 size_t size()
const;
18677 void set_rand_engine(std::shared_ptr< std::mt19937 > & eng);
18678 std::shared_ptr< std::mt19937 > & get_rand_endgine();
18679 void seed(
size_t s);
18680 void set_rand_norm(epiworld_double mean, epiworld_double sd);
18681 void set_rand_unif(epiworld_double a, epiworld_double b);
18682 void set_rand_exp(epiworld_double lambda);
18683 void set_rand_gamma(epiworld_double alpha, epiworld_double beta);
18684 void set_rand_lognormal(epiworld_double mean, epiworld_double shape);
18685 void set_rand_binom(
int n, epiworld_double p);
18686 void set_rand_nbinom(
int n, epiworld_double p);
18687 void set_rand_geom(epiworld_double p);
18688 void set_rand_poiss(epiworld_double lambda);
18689 epiworld_double runif();
18690 epiworld_double runif(epiworld_double a, epiworld_double b);
18691 epiworld_double rnorm();
18692 epiworld_double rnorm(epiworld_double mean, epiworld_double sd);
18693 epiworld_double rgamma();
18694 epiworld_double rgamma(epiworld_double alpha, epiworld_double beta);
18695 epiworld_double rexp();
18696 epiworld_double rexp(epiworld_double lambda);
18697 epiworld_double rlognormal();
18698 epiworld_double rlognormal(epiworld_double mean, epiworld_double shape);
18700 int rbinom(
int n, epiworld_double p);
18702 int rnbinom(
int n, epiworld_double p);
18704 int rgeom(epiworld_double p);
18706 int rpoiss(epiworld_double lambda);
18724 void rm_virus(
size_t virus_pos);
18725 void rm_tool(
size_t tool_pos);
18726 void rm_entity(
size_t entity_id);
18739 void load_agents_entities_ties(std::string fn,
int skip);
18744 void load_agents_entities_ties(
18745 const std::vector<int> & agents_ids,
18746 const std::vector<int> & entities_ids
18749 void load_agents_entities_ties(
18750 const int * agents_id,
18751 const int * entities_id,
18765 void agents_from_adjlist(
18769 bool directed =
false
18772 void agents_from_edgelist(
18773 const std::vector< int > & source,
18774 const std::vector< int > & target,
18779 void agents_from_adjlist(
AdjList al);
18781 bool is_directed()
const;
18783 std::vector< Agent<TSeq> > & get_agents();
18787 std::vector< epiworld_fast_uint > get_agents_states()
const;
18789 std::vector< Viruses_const<TSeq> > get_agents_viruses()
const;
18791 std::vector< Viruses<TSeq> > get_agents_viruses();
18793 std::vector< Entity<TSeq> > & get_entities();
18795 Entity<TSeq> & get_entity(
size_t entity_id,
int * entity_pos =
nullptr);
18798 epiworld_fast_uint n = 1000,
18799 epiworld_fast_uint k = 5,
18801 epiworld_double p = .01
18803 void agents_empty_graph(epiworld_fast_uint n = 1000);
18816 void update_state();
18817 void mutate_virus();
18820 epiworld_fast_uint ndays,
18824 epiworld_fast_uint ndays,
18825 epiworld_fast_uint nexperiments,
18827 std::function<
void(
size_t,
Model<TSeq>*)> fun = make_save_run<TSeq>(),
18829 bool verbose =
true,
18834 size_t get_n_viruses()
const;
18835 size_t get_n_tools()
const;
18836 epiworld_fast_uint get_ndays()
const;
18837 epiworld_fast_uint get_n_replicates()
const;
18838 size_t get_n_entities()
const;
18839 void set_ndays(epiworld_fast_uint ndays);
18840 bool get_verbose()
const;
18858 void set_rewire_prop(epiworld_double prop);
18859 epiworld_double get_rewire_prop()
const;
18876 std::string fn_virus_info,
18877 std::string fn_virus_hist,
18878 std::string fn_tool_info,
18879 std::string fn_tool_hist,
18880 std::string fn_total_hist,
18881 std::string fn_transmission,
18882 std::string fn_transition,
18883 std::string fn_reproductive_number,
18884 std::string fn_generation_time
18898 void write_edgelist(
18902 void write_edgelist(
18903 std::vector< int > & source,
18904 std::vector< int > & target
18908 std::map<std::string, epiworld_double> & params();
18921 virtual void reset();
18922 const Model<TSeq> & print(
bool lite =
false)
const;
18939 void add_state(std::string lab, UpdateFun<TSeq> fun =
nullptr);
18940 const std::vector< std::string > & get_states()
const;
18941 size_t get_n_states()
const;
18942 const std::vector< UpdateFun<TSeq> > & get_state_fun()
const;
18943 void print_state_codes()
const;
18955 std::vector< double > ,
18994 epiworld_double add_param(
18995 epiworld_double initial_val, std::string pname,
bool overwrite =
false
18997 Model<TSeq> & read_params(std::string fn,
bool overwrite =
false);
18998 epiworld_double get_param(epiworld_fast_uint k);
18999 epiworld_double get_param(std::string pname);
19001 void set_param(std::string pname, epiworld_double val);
19003 epiworld_double par(std::string pname)
const;
19007 std::string unit =
"auto",
19008 epiworld_double * last_elapsed =
nullptr,
19009 epiworld_double * total_elapsed =
nullptr,
19010 std::string * unit_abbr =
nullptr,
19020 void set_user_data(std::vector< std::string > names);
19021 void add_user_data(epiworld_fast_uint j, epiworld_double x);
19022 void add_user_data(std::vector< epiworld_double > x);
19037 void add_globalevent(
19039 std::string name =
"A global action",
19043 void add_globalevent(
19050 void rm_globalevent(std::string name);
19051 void rm_globalevent(
size_t i);
19053 void run_globalevents();
19055 void clear_state_set();
19067 bool is_queuing_on()
const;
19078 void set_susceptibility_reduction_mixer(MixerFun<TSeq> fun);
19079 void set_transmission_reduction_mixer(MixerFun<TSeq> fun);
19080 void set_recovery_enhancer_mixer(MixerFun<TSeq> fun);
19081 void set_death_reduction_mixer(MixerFun<TSeq> fun);
19084 const std::vector< VirusPtr<TSeq> > & get_viruses()
const;
19085 const std::vector< ToolPtr<TSeq> > & get_tools()
const;
19100 void set_agents_data(
double * data_,
size_t ncols_);
19101 double * get_agents_data();
19102 size_t get_agents_data_ncols()
const;
19109 void set_name(std::string name);
19110 std::string get_name()
const;
19112 bool operator==(
const Model<TSeq> & other)
const;
19113 bool operator!=(
const Model<TSeq> & other)
const {
return !operator==(other);};
19130 DiagramType diagram_type = DiagramType::Mermaid,
19131 const std::string & fn_output =
"",
19210 m->array_double_tmp[nviruses_tmp] =
19211 (1.0 - p->get_susceptibility_reduction(v, m)) *
19212 v->get_prob_infecting(m) *
19213 (1.0 - neighbor->get_transmission_reduction(v, m))
19216 m->array_virus_tmp[nviruses_tmp++] = &(*v);
19221 if (nviruses_tmp == 0u)
19225 int which = roulette(nviruses_tmp, m);
19230 p->set_virus(*m->array_virus_tmp[which], m);
19240 std::shared_ptr<std::vector<bool>> exclude_agent_bool =
19241 std::make_shared<std::vector<bool>>(0);
19243 std::shared_ptr<std::vector<epiworld_fast_uint>> exclude_agent_bool_idx =
19244 std::make_shared<std::vector<epiworld_fast_uint>>(exclude);
19251 if (exclude_agent_bool->size() == 0u)
19254 exclude_agent_bool->resize(m->get_states().size(),
false);
19255 for (
auto s : *exclude_agent_bool_idx)
19257 if (s >= exclude_agent_bool->size())
19258 throw std::logic_error(
19259 std::string(
"You are trying to exclude a state that is out of range: ") +
19260 std::to_string(s) + std::string(
". There are only ") +
19261 std::to_string(exclude_agent_bool->size()) +
19262 std::string(
" states in the model.")
19265 exclude_agent_bool->operator[](s) =
true;
19271 if (p->get_virus() !=
nullptr)
19272 throw std::logic_error(
19273 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense! ") +
19274 std::string(
"Agent id ") + std::to_string(p->get_id()) +
19275 std::string(
" has a virus.")
19279 size_t nviruses_tmp = 0u;
19280 for (
auto & neighbor: p->get_neighbors())
19284 if (exclude_agent_bool->operator[](neighbor->get_state()))
19287 auto & v = neighbor->get_virus();
19293 m->array_double_tmp[nviruses_tmp] =
19294 (1.0 - p->get_susceptibility_reduction(v, m)) *
19295 v->get_prob_infecting(m) *
19296 (1.0 - neighbor->get_transmission_reduction(v, m))
19299 m->array_virus_tmp[nviruses_tmp++] = &(*v);
19304 if (nviruses_tmp == 0u)
19308 int which = roulette(nviruses_tmp, m);
19313 p->set_virus(*m->array_virus_tmp[which], m);
19338 template<
typename TSeq = EPI_DEFAULT_TSEQ>
19340 std::vector< epiworld_fast_uint > exclude = {}
19343 if (exclude.size() == 0u)
19349 if (p->get_virus() !=
nullptr)
19350 throw std::logic_error(
19351 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense! ") +
19352 std::string(
"Agent id ") + std::to_string(p->get_id()) +
19353 std::string(
" has a virus.")
19357 size_t nviruses_tmp = 0u;
19358 for (
auto & neighbor: p->get_neighbors())
19361 if (neighbor->get_virus() ==
nullptr)
19364 auto & v = neighbor->get_virus();
19367 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
19368 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
19372 m->array_double_tmp[nviruses_tmp] =
19373 (1.0 - p->get_susceptibility_reduction(v, m)) *
19374 v->get_prob_infecting(m) *
19375 (1.0 - neighbor->get_transmission_reduction(v, m))
19378 m->array_virus_tmp[nviruses_tmp++] = &(*v);
19383 if (nviruses_tmp == 0u)
19387 int which = roulette(nviruses_tmp, m);
19392 return m->array_virus_tmp[which];
19402 std::shared_ptr<std::vector<bool>> exclude_agent_bool =
19403 std::make_shared<std::vector<bool>>(0);
19405 std::shared_ptr<std::vector<epiworld_fast_uint>> exclude_agent_bool_idx =
19406 std::make_shared<std::vector<epiworld_fast_uint>>(exclude);
19413 if (exclude_agent_bool->size() == 0u)
19416 exclude_agent_bool->resize(m->get_states().size(),
false);
19417 for (
auto s : *exclude_agent_bool_idx)
19419 if (s >= exclude_agent_bool->size())
19420 throw std::logic_error(
19421 std::string(
"You are trying to exclude a state that is out of range: ") +
19422 std::to_string(s) + std::string(
". There are only ") +
19423 std::to_string(exclude_agent_bool->size()) +
19424 std::string(
" states in the model.")
19427 exclude_agent_bool->operator[](s) =
true;
19433 if (p->get_virus() !=
nullptr)
19434 throw std::logic_error(
19435 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense! ") +
19436 std::string(
"Agent id ") + std::to_string(p->get_id()) +
19437 std::string(
" has a virus.")
19441 size_t nviruses_tmp = 0u;
19442 for (
auto & neighbor: p->get_neighbors())
19446 if (exclude_agent_bool->operator[](neighbor->get_state()))
19449 if (neighbor->get_virus() ==
nullptr)
19452 auto & v = neighbor->get_virus();
19455 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
19456 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
19460 m->array_double_tmp[nviruses_tmp] =
19461 (1.0 - p->get_susceptibility_reduction(v, m)) *
19462 v->get_prob_infecting(m) *
19463 (1.0 - neighbor->get_transmission_reduction(v, m))
19466 m->array_virus_tmp[nviruses_tmp++] = &(*v);
19471 if (nviruses_tmp == 0u)
19475 int which = roulette(nviruses_tmp, m);
19480 return m->array_virus_tmp[which];
19506 template<
typename TSeq = EPI_DEFAULT_TSEQ>
19510 if (p->get_virus() !=
nullptr)
19511 throw std::logic_error(
19512 std::string(
"Using the -default_update_susceptible- on agents WITH viruses makes no sense!") +
19513 std::string(
"Agent id ") + std::to_string(p->
get_id()) +
19514 std::string(
" has a virus.")
19518 size_t nviruses_tmp = 0u;
19519 for (
auto & neighbor: p->get_neighbors())
19522 int _vcount_neigh = 0;
19525 if (neighbor->get_virus() ==
nullptr)
19528 auto & v = neighbor->get_virus();
19531 if (nviruses_tmp >= m->array_virus_tmp.size())
19532 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
19536 m->array_double_tmp[nviruses_tmp] =
19537 (1.0 - p->get_susceptibility_reduction(v, m)) *
19538 v->get_prob_infecting(m) *
19539 (1.0 - neighbor->get_transmission_reduction(v, m))
19542 m->array_virus_tmp[nviruses_tmp++] = &(*v);
19546 (m->array_double_tmp[nviruses_tmp - 1] < 0.0) |
19547 (m->array_double_tmp[nviruses_tmp - 1] > 1.0)
19551 "[epi-debug] Agent %i's virus %i has transmission prob outside of [0, 1]: %.4f!\n",
19552 static_cast<int>(neighbor->get_id()),
19553 static_cast<int>(_vcount_neigh++),
19554 m->array_double_tmp[nviruses_tmp - 1]
19563 if (nviruses_tmp == 0u)
19567 m->get_db().n_transmissions_potential++;
19571 int which = roulette(nviruses_tmp, m);
19577 m->get_db().n_transmissions_today++;
19580 return m->array_virus_tmp[which];
21932 template<
typename TSeq>
21934 epiworld_fast_uint tool_idx,
21936 epiworld_fast_int state_new,
21937 epiworld_fast_int queue
21941 if (tool_idx >= n_tools)
21942 throw std::range_error(
21943 "The Tool you want to remove is out of range. This Agent only has " +
21944 std::to_string(n_tools) +
" tools."
21948 this,
nullptr, tools[tool_idx],
nullptr, state_new, queue, default_rm_tool<TSeq>, -1, -1
21953 template<
typename TSeq>
21955 ToolPtr<TSeq> & tool,
21957 epiworld_fast_int state_new,
21958 epiworld_fast_int queue
21962 if (tool->agent !=
this)
21963 throw std::logic_error(
"Cannot remove a virus from another agent!");
21966 this,
nullptr, tool,
nullptr, state_new, queue, default_rm_tool<TSeq>, -1, -1
21971 template<
typename TSeq>
21974 epiworld_fast_int state_new,
21975 epiworld_fast_int queue
21979 if (virus ==
nullptr)
21980 throw std::logic_error(
21981 "There is no virus to remove here!"
21984 if (state_new == -99)
21985 virus->get_state(
nullptr, &state_new,
nullptr);
21988 virus->get_queue(
nullptr, &queue,
nullptr);
21991 this, virus,
nullptr,
nullptr,
21994 default_rm_virus<TSeq>, -1, -1
21999 template<
typename TSeq>
22001 epiworld_fast_uint entity_idx,
22003 epiworld_fast_int state_new,
22004 epiworld_fast_int queue
22008 if (entity_idx >= n_entities)
22009 throw std::range_error(
22010 "The Entity you want to remove is out of range. This Agent only has " +
22011 std::to_string(n_entities) +
" entitites."
22013 else if (n_entities == 0u)
22014 throw std::logic_error(
22015 "There is entity to remove here!"
22022 &model->get_entity(entity_idx),
22025 default_rm_entity<TSeq>,
22026 entities_locations[entity_idx],
22031 template<
typename TSeq>
22035 epiworld_fast_int state_new,
22036 epiworld_fast_int queue
22041 int entity_idx = -1;
22042 for (
size_t i = 0u; i < n_entities; ++i)
22044 if (
static_cast<int>(entities[i]) == entity.get_id())
22051 if (entity_idx == -1)
22052 throw std::logic_error(
22053 std::string(
"The agent ") +
22054 std::to_string(
id) +
22055 std::string(
" is not associated with entity \"") +
22056 entity.get_name() +
22064 &model->entities[entity.get_id()],
22067 default_rm_entity<TSeq>,
22068 entities_locations[entity_idx],
22073 template<
typename TSeq>
22080 this, virus,
nullptr,
nullptr,
22081 virus->state_removed,
22082 virus->queue_removed,
22083 default_rm_virus<TSeq>, -1, -1
22088 template<
typename TSeq>
22094 return model->susceptibility_reduction_mixer(
this, v, model);
22097 template<
typename TSeq>
22102 return model->transmission_reduction_mixer(
this, v, model);
22105 template<
typename TSeq>
22110 return model->recovery_enhancer_mixer(
this, v, model);
22113 template<
typename TSeq>
22118 return model->death_reduction_mixer(
this, v, model);
22121 template<
typename TSeq>
22127 template<
typename TSeq>
22132 template<
typename TSeq>
22138 template<
typename TSeq>
22143 template<
typename TSeq>
22148 template<
typename TSeq>
22151 return tools.at(i);
22154 template<
typename TSeq>
22160 template<
typename TSeq>
22168 template<
typename TSeq>
22175 bool found =
false;
22177 if (neighbors ==
nullptr)
22179 neighbors =
new std::vector< size_t >();
22180 neighbors_locations =
new std::vector< size_t >();
22183 if (check_source && neighbors)
22186 for (
auto & n: *neighbors)
22187 if (
static_cast<int>(n) == p.
get_id())
22202 neighbors_locations->push_back(p.get_n_neighbors());
22203 neighbors->push_back(p.
get_id());
22210 if (check_target && p.neighbors)
22213 for (
auto & n: *p.neighbors)
22214 if (
static_cast<int>(n) ==
id)
22225 if (p.neighbors ==
nullptr)
22227 p.neighbors =
new std::vector< size_t >();
22228 p.neighbors_locations =
new std::vector< size_t >();
22231 p.neighbors_locations->push_back(n_neighbors - 1);
22232 p.neighbors->push_back(
id);
22240 template<
typename TSeq>
22248 if (n_this >= n_neighbors)
22249 throw std::range_error(
22250 "The neighbor you want to swap is out of range. This Agent only has " +
22251 std::to_string(n_neighbors) +
" neighbors."
22253 if (n_other >= other.n_neighbors)
22254 throw std::range_error(
22255 "The neighbor you want to swap is out of range. This Agent only has " +
22256 std::to_string(other.n_neighbors) +
" neighbors."
22260 auto & pop = model->population;
22261 auto & neigh_this = pop[(*neighbors)[n_this]];
22262 auto & neigh_other = pop[(*other.neighbors)[n_other]];
22265 size_t loc_this_in_neigh = (*neighbors_locations)[n_this];
22266 size_t loc_other_in_neigh = (*other.neighbors_locations)[n_other];
22269 std::swap((*neighbors)[n_this], (*other.neighbors)[n_other]);
22271 if (!model->directed)
22274 (*neigh_this.neighbors)[loc_this_in_neigh],
22275 (*neigh_other.neighbors)[loc_other_in_neigh]
22279 std::swap((*neighbors_locations)[n_this], (*other.neighbors_locations)[n_other]);
22282 (*neigh_this.neighbors_locations)[loc_this_in_neigh],
22283 (*neigh_other.neighbors_locations)[loc_other_in_neigh]
22289 template<
typename TSeq>
22292 std::vector< Agent<TSeq> * > res(n_neighbors,
nullptr);
22293 for (
size_t i = 0u; i < n_neighbors; ++i)
22294 res[i] = &model->population[(*neighbors)[i]];
22299 template<
typename TSeq>
22302 return n_neighbors;
22305 template<
typename TSeq>
22308 epiworld_fast_uint new_state,
22309 epiworld_fast_int queue
22314 this,
nullptr,
nullptr,
nullptr, new_state, queue,
22315 default_change_state<TSeq>, -1, -1
22322 template<
typename TSeq>
22327 template<
typename TSeq>
22331 this->virus =
nullptr;
22333 this->tools.clear();
22336 this->entities.clear();
22337 this->entities_locations.clear();
22338 this->n_entities = 0u;
22341 this->state_prev = 0u;
22343 this->state_last_changed = -1;
22347 template<
typename TSeq>
22351 for (
auto & tool : tools)
22352 if (tool->get_id() ==
static_cast<int>(t))
22359 template<
typename TSeq>
22363 for (
auto & tool : tools)
22364 if (tool->get_name() == name)
22371 template<
typename TSeq>
22375 return has_tool(tool.get_id());
22379 template<
typename TSeq>
22382 if (virus->get_id() ==
static_cast<int>(t))
22388 template<
typename TSeq>
22392 if (virus->get_name() == name)
22399 template<
typename TSeq>
22403 return has_virus(virus.get_id());
22407 template<
typename TSeq>
22411 for (
auto & entity : entities)
22419 template<
typename TSeq>
22423 for (
auto & entity : entities)
22424 if (model->get_entity(entity).get_name() == name)
22431 template<
typename TSeq>
22440 "Agent: %i, state: %s (%i), Has virus: %s, NTools: %ii NNeigh: %i\n",
22441 static_cast<int>(
id),
22442 model->states_labels[state].c_str(),
22443 static_cast<int>(state),
22444 virus ==
nullptr ? std::string(
"no").c_str() : std::string(
"yes").c_str(),
22445 static_cast<int>(n_tools),
22446 static_cast<int>(n_neighbors)
22451 printf_epiworld(
"Information about agent id %i\n",
22452 static_cast<int>(this->
id));
22453 printf_epiworld(
" State : %s (%i)\n",
22454 model->states_labels[state].c_str(),
static_cast<int>(state));
22455 printf_epiworld(
" Has virus : %s\n", virus ==
nullptr ?
22456 std::string(
"no").c_str() : std::string(
"yes").c_str());
22457 printf_epiworld(
" Tool count : %i\n",
static_cast<int>(n_tools));
22458 printf_epiworld(
" Neigh. count : %i\n",
static_cast<int>(n_neighbors));
22460 size_t nfeats = model->get_agents_data_ncols();
22465 "This model includes features (%i): [ ",
22466 static_cast<int>(nfeats)
22469 int max_to_show =
static_cast<int>((nfeats > 10)? 10 : nfeats);
22471 for (
int k = 0; k < max_to_show; ++k)
22473 printf_epiworld(
"%.2f", this->
operator[](k));
22475 if (k != (max_to_show - 1))
22477 printf_epiworld(
", ");
22479 printf_epiworld(
" ]\n");
22492 template<
typename TSeq>
22496 if (model->agents_data_ncols <= j)
22497 throw std::logic_error(
"The requested feature of the agent is out of range.");
22499 return *(model->agents_data + j * model->size() + id);
22503 template<
typename TSeq>
22506 return *(model->agents_data + j * model->size() + id);
22509 template<
typename TSeq>
22513 if (model->agents_data_ncols <= j)
22514 throw std::logic_error(
"The requested feature of the agent is out of range.");
22516 return *(model->agents_data + j * model->size() + id);
22520 template<
typename TSeq>
22523 return *(model->agents_data + j * model->size() + id);
22526 template<
typename TSeq>
22532 template<
typename TSeq>
22538 template<
typename TSeq>
22541 if (n_entities == 0)
22542 throw std::range_error(
"Agent id " + std::to_string(
id) +
" has no entities.");
22544 if (i >= n_entities)
22545 throw std::range_error(
"Trying to get to an agent's entity outside of the range.");
22547 return model->get_entity(entities[i]);
22550 template<
typename TSeq>
22553 if (n_entities == 0)
22554 throw std::range_error(
"Agent id " + std::to_string(
id) +
" has no entities.");
22556 if (i >= n_entities)
22557 throw std::range_error(
"Trying to get to an agent's entity outside of the range.");
22559 return model->get_entity(entities[i]);
22562 template<
typename TSeq>
22568 template<
typename TSeq>
22572 EPI_DEBUG_FAIL_AT_TRUE(
22573 n_neighbors != other.n_neighbors,
22574 "Agent:: n_eighbors don't match"
22578 for (
size_t i = 0u; i < n_neighbors; ++i)
22580 EPI_DEBUG_FAIL_AT_TRUE(
22581 (*neighbors)[i] != (*other.neighbors)[i],
22582 "Agent:: neighbor[i] don't match"
22586 EPI_DEBUG_FAIL_AT_TRUE(
22587 n_entities != other.n_entities,
22588 "Agent:: n_entities don't match"
22592 for (
size_t i = 0u; i < n_entities; ++i)
22594 EPI_DEBUG_FAIL_AT_TRUE(
22595 entities[i] != other.entities[i],
22596 "Agent:: entities[i] don't match"
22600 EPI_DEBUG_FAIL_AT_TRUE(
22601 state != other.state,
22602 "Agent:: state don't match"
22606 EPI_DEBUG_FAIL_AT_TRUE(
22607 state_prev != other.state_prev,
22608 "
Agent:: state_prev don't match"
22617 EPI_DEBUG_FAIL_AT_TRUE(
22618 ((virus ==
nullptr) && (other.virus !=
nullptr)) ||
22619 ((virus !=
nullptr) && (other.virus ==
nullptr)),
22620 "
Agent:: virus don't match"
22623 if ((virus !=
nullptr) && (other.virus !=
nullptr))
22625 EPI_DEBUG_FAIL_AT_TRUE(
22626 *virus != *other.virus,
22627 "Agent:: virus doesn't match"
22631 EPI_DEBUG_FAIL_AT_TRUE(n_tools != other.n_tools,
"Agent:: n_tools don't match")
22633 for (
size_t i = 0u; i < n_tools; ++i)
22636 EPI_DEBUG_FAIL_AT_TRUE(
22637 tools[i] != other.tools[i],
22638 "Agent:: tools[i] don't match"
24094 create_init_function_sir<TSeq>(proportions_)
24272 create_init_function_seir<TSeq>(proportions_)
24448 epiworld_double tmp_transmission =
24449 (1.0 - p->get_susceptibility_reduction(v, m)) *
24450 v->get_prob_infecting(m) *
24451 (1.0 - neighbor->get_transmission_reduction(v, m))
24454 m->array_double_tmp[nviruses_tmp] = tmp_transmission;
24455 m->array_virus_tmp[nviruses_tmp++] = &(*v);
24459 if (nviruses_tmp == 0)
24463 int which = roulette(nviruses_tmp, m);
24468 p->set_virus(*m->array_virus_tmp[which], m);
24474 epiworld::UpdateFun<TSeq> surveillance_update_exposed =
24481 epiworld::VirusPtr<TSeq> & v = p->get_virus();
24482 epiworld_double p_die = v->get_prob_death(m) * (1.0 - p->get_death_reduction(v, m));
24484 epiworld_fast_uint days_since_exposed = m->
today() - v->get_date();
24485 epiworld_fast_uint state = p->get_state();
24489 if (dat[p->
get_id()] < 0)
24491 epiworld_double latent_days = m->rgamma(
24492 m->par(
"Latent period"), 1.0
24495 dat[p->
get_id() * 2u] = latent_days;
24497 dat[p->
get_id() * 2u + 1u] =
24498 m->rgamma(m->par(
"Infect period"), 1.0) +
24503 if (days_since_exposed <= dat[p->
get_id() * 2u])
24507 if (days_since_exposed >= dat[p->
get_id() * 2u + 1u])
24518 if (EPI_RUNIF() < m->par(
"Prob of symptoms"))
24528 if (EPI_RUNIF() < p_die)
24538 std::vector< unsigned int > exposed_state = {
24540 SYMPTOMATIC_ISOLATED,
24542 ASYMPTOMATIC_ISOLATED,
24546 epiworld::GlobalFun<TSeq> surveillance_program =
24553 std::binomial_distribution<> bdist(m->size(), m->par(
"Surveilance prob."));
24554 int nsampled = bdist(*m->get_rand_endgine());
24556 int to_go = nsampled + 1;
24558 epiworld_double ndetected = 0.0;
24559 epiworld_double ndetected_asympt = 0.0;
24562 std::vector< bool > sampled(m->size(),
false);
24564 while (to_go-- > 0)
24568 epiworld_fast_uint i =
static_cast<epiworld_fast_uint
>(std::floor(EPI_RUNIF() * m->size()));
24577 if (epiworld::IN(p->get_state(), exposed_state ))
24583 ndetected_asympt += 1.0;
24596 std::vector< int > totals;
24597 m->get_db().get_today_total(
nullptr,&totals);
24600 static_cast<epiworld_double
>(nsampled),
24610 model.add_state(
"Susceptible", surveillance_update_susceptible);
24611 model.add_state(
"Latent", surveillance_update_exposed);
24612 model.add_state(
"Symptomatic", surveillance_update_exposed);
24613 model.add_state(
"Symptomatic isolated", surveillance_update_exposed);
24614 model.add_state(
"Asymptomatic", surveillance_update_exposed);
24615 model.add_state(
"Asymptomatic isolated", surveillance_update_exposed);
24616 model.add_state(
"Recovered");
24617 model.add_state(
"Removed");
24620 model.add_param(latent_period,
"Latent period");
24621 model.add_param(infect_period,
"Infect period");
24622 model.add_param(prob_symptoms,
"Prob of symptoms");
24623 model.add_param(surveillance_prob,
"Surveilance prob.");
24624 model.add_param(efficacy_vax,
"Vax efficacy");
24625 model.add_param(prop_vax_redux_transm,
"Vax redux transmission");
24626 model.add_param(prob_transmission,
"Prob of transmission");
24627 model.add_param(prob_death,
"Prob. death");
24628 model.add_param(prob_noreinfect,
"Prob. no reinfect");
24632 covid.set_state(LATENT, RECOVERED, REMOVED);
24633 covid.set_post_immunity(&model(
"Prob. no reinfect"));
24634 covid.set_prob_death(&model(
"Prob. death"));
24636 epiworld::VirusFun<TSeq> ptransmitfun = [](
24640 ) -> epiworld_double
24643 epiworld_fast_uint s = p->get_state();
24645 return static_cast<epiworld_double
>(0.0);
24647 return static_cast<epiworld_double
>(0.0);
24649 return static_cast<epiworld_double
>(0.0);
24652 return m->par(
"Prob of transmission");
24655 covid.set_prob_infecting_fun(ptransmitfun);
24657 model.add_virus(covid);
24659 model.set_user_data({
"nsampled",
"ndetected",
"ndetected_asympt",
"nasymptomatic"});
24660 model.add_globalevent(surveillance_program,
"Surveilance program", -1);
24664 vax.set_susceptibility_reduction(&model(
"Vax efficacy"));
24665 vax.set_transmission_reduction(&model(
"Vax redux transmission"));
24667 model.add_tool(vax);
24669 model.set_name(
"Surveillance");
24675 template<
typename TSeq>
24677 const std::string & vname,
24678 epiworld_fast_uint prevalence,
24679 epiworld_double efficacy_vax,
24680 epiworld_double latent_period,
24681 epiworld_double infect_period,
24682 epiworld_double prob_symptoms,
24683 epiworld_double prop_vaccinated,
24684 epiworld_double prop_vax_redux_transm,
24685 epiworld_double prop_vax_redux_infect,
24686 epiworld_double surveillance_prob,
24687 epiworld_double prob_transmission,
24688 epiworld_double prob_death,
24689 epiworld_double prob_noreinfect
24702 prop_vax_redux_transm,
24703 prop_vax_redux_infect,
24806 std::vector< double > generation_time_expected(
24807 int max_days = 200,
24808 int max_contacts = 200
24813 template<
typename TSeq>
24818 infected.reserve(this->size());
24820 for (
auto & p : this->get_agents())
24824 infected.push_back(&p);
24829 this->get_n_infected(),
24838 template<
typename TSeq>
24840 epiworld_fast_uint ndays,
24850 template<
typename TSeq>
24856 this->update_infected();
24862 template<
typename TSeq>
24884 template<
typename TSeq>
24887 const std::string & vname,
24888 epiworld_fast_uint n,
24889 epiworld_double prevalence,
24890 epiworld_double contact_rate,
24891 epiworld_double transmission_rate,
24892 epiworld_double recovery_rate
24897 epiworld::UpdateFun<TSeq> update_susceptible = [](
24902 int ndraw = m->rbinom();
24908 size_t ninfected = model->get_n_infected();
24911 int nviruses_tmp = 0;
24912 for (
int i = 0; i < ndraw; ++i)
24915 int which =
static_cast<int>(
24916 std::floor(ninfected * m->runif())
24926 if (which ==
static_cast<int>(ninfected))
24936 if (neighbor.get_virus() ==
nullptr)
24939 auto & v = neighbor.get_virus();
24942 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
24943 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
24947 m->array_double_tmp[nviruses_tmp] =
24948 (1.0 - p->get_susceptibility_reduction(v, m)) *
24949 v->get_prob_infecting(m) *
24950 (1.0 - neighbor.get_transmission_reduction(v, m))
24953 m->array_virus_tmp[nviruses_tmp++] = &(*v);
24958 if (nviruses_tmp == 0u)
24962 int which = roulette(nviruses_tmp, m);
24967 p->set_virus(*m->array_virus_tmp[which], m);
24974 epiworld::UpdateFun<TSeq> update_infected = [](
24978 auto state = p->get_state();
24985 epiworld_fast_uint n_events = 0u;
24987 m->array_double_tmp[n_events++] =
24988 1.0 - (1.0 - p->get_virus()->get_prob_recovery(m)) *
24989 (1.0 - p->get_recovery_enhancer(p->get_virus(), m));
24992 if (n_events == 0u)
24995 "[epi-debug] agent %i has 0 possible events!!\n",
24996 static_cast<int>(p->
get_id())
24998 throw std::logic_error(
"Zero events in exposed.");
25001 if (n_events == 0u)
25007 int which = roulette(n_events, m);
25018 throw std::logic_error(
25019 "This function can only be applied to infected individuals. (SIR)"
25027 model.add_state(
"Susceptible", update_susceptible);
25028 model.add_state(
"Infected", update_infected);
25029 model.add_state(
"Recovered");
25032 model.add_param(contact_rate,
"Contact rate");
25033 model.add_param(transmission_rate,
"Transmission rate");
25034 model.add_param(recovery_rate,
"Recovery rate");
25041 model->update_infected();
25046 model.add_globalevent(update,
"Update infected individuals");
25050 virus.set_state(1, 2, 2);
25051 virus.set_prob_infecting(&model(
"Transmission rate"));
25052 virus.set_prob_recovery(&model(
"Recovery rate"));
25054 model.add_virus(virus);
25056 model.queuing_off();
25058 model.agents_empty_graph(n);
25060 model.set_name(
"Susceptible-Infected-Removed (SIR) (connected)");
25066 template<
typename TSeq>
25068 const std::string & vname,
25069 epiworld_fast_uint n,
25070 epiworld_double prevalence,
25071 epiworld_double contact_rate,
25072 epiworld_double transmission_rate,
25073 epiworld_double recovery_rate
25091 template<
typename TSeq>
25093 std::vector< double > proportions_,
25098 create_init_function_sir<TSeq>(proportions_)
25105 template<
typename TSeq>
25113 std::vector< int > h_date;
25114 std::vector< std::string > h_state;
25115 std::vector< int > h_counts;
25117 this_const->get_db().get_hist_total(
25124 std::vector< double > S(this_const->get_ndays(), 0.0);
25125 for (
size_t i = 0; i < h_date.size(); ++i)
25127 if (h_state[i] ==
"Susceptible")
25128 S[h_date[i]] += h_counts[i];
25134 std::vector< double > gen_times(this_const->get_ndays(), 1.0);
25135 double p_c = this_const->par(
"Contact rate")/this_const->size();
25136 double p_i = this_const->par(
"Transmission rate");
25137 double p_r = this_const->par(
"Recovery rate");
25138 for (
size_t i = 0u; i < this_const->get_ndays(); ++i)
25140 gen_times[i] = gen_int_mean(
25241 std::vector< double > generation_time_expected(
25242 int max_days = 200,
25243 int max_contacts = 200
25248 template<
typename TSeq>
25253 infected.reserve(this->size());
25255 for (
auto & p : this->get_agents())
25259 infected.push_back(&p);
25264 this->get_n_infected(),
25273 template<
typename TSeq>
25275 epiworld_fast_uint ndays,
25286 template<
typename TSeq>
25291 this->update_infected();
25297 template<
typename TSeq>
25319 template<
typename TSeq>
25322 const std::string & vname,
25323 epiworld_fast_uint n,
25324 epiworld_double prevalence,
25325 epiworld_double contact_rate,
25326 epiworld_double transmission_rate,
25327 epiworld_double avg_incubation_days,
25328 epiworld_double recovery_rate
25333 epiworld::UpdateFun<TSeq> update_susceptible = [](
25339 int ndraw = m->rbinom();
25345 size_t ninfected = model->get_n_infected();
25348 int nviruses_tmp = 0;
25349 for (
int i = 0; i < ndraw; ++i)
25352 int which =
static_cast<int>(
25353 std::floor(ninfected * m->runif())
25363 if (which ==
static_cast<int>(ninfected))
25373 auto & v = neighbor.get_virus();
25376 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
25377 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
25381 m->array_double_tmp[nviruses_tmp] =
25382 (1.0 - p->get_susceptibility_reduction(v, m)) *
25383 v->get_prob_infecting(m) *
25384 (1.0 - neighbor.get_transmission_reduction(v, m))
25387 m->array_virus_tmp[nviruses_tmp++] = &(*v);
25392 if (nviruses_tmp == 0u)
25396 int which = roulette(nviruses_tmp, m);
25402 *m->array_virus_tmp[which],
25411 epiworld::UpdateFun<TSeq> update_infected = [](
25415 auto state = p->get_state();
25421 auto & v = p->get_virus();
25424 if (m->runif() < 1.0/(v->get_incubation(m)))
25438 epiworld_fast_uint n_events = 0u;
25439 const auto & v = p->get_virus();
25442 m->array_double_tmp[n_events++] =
25443 1.0 - (1.0 - v->get_prob_recovery(m)) * (1.0 - p->get_recovery_enhancer(v, m));
25446 if (n_events == 0u)
25449 "[epi-debug] agent %i has 0 possible events!!\n",
25450 static_cast<int>(p->
get_id())
25452 throw std::logic_error(
"Zero events in exposed.");
25455 if (n_events == 0u)
25461 int which = roulette(n_events, m);
25472 throw std::logic_error(
"This function can only be applied to exposed or infected individuals. (SEIR)") ;
25479 model.add_param(contact_rate,
"Contact rate");
25480 model.add_param(transmission_rate,
"Prob. Transmission");
25481 model.add_param(recovery_rate,
"Prob. Recovery");
25482 model.add_param(avg_incubation_days,
"Avg. Incubation days");
25485 model.add_state(
"Susceptible", update_susceptible);
25486 model.add_state(
"Exposed", update_infected);
25487 model.add_state(
"Infected", update_infected);
25488 model.add_state(
"Recovered");
25491 epiworld::GlobalFun<TSeq> update = [](
25498 model->update_infected();
25504 model.add_globalevent(update,
"Update infected individuals");
25515 virus.set_prob_infecting(&model(
"Prob. Transmission"));
25516 virus.set_prob_recovery(&model(
"Prob. Recovery"));
25517 virus.set_incubation(&model(
"Avg. Incubation days"));
25519 model.add_virus(virus);
25521 model.queuing_off();
25524 model.agents_empty_graph(n);
25526 model.set_name(
"Susceptible-Exposed-Infected-Removed (SEIR) (connected)");
25532 template<
typename TSeq>
25534 const std::string & vname,
25535 epiworld_fast_uint n,
25536 epiworld_double prevalence,
25537 epiworld_double contact_rate,
25538 epiworld_double transmission_rate,
25539 epiworld_double avg_incubation_days,
25540 epiworld_double recovery_rate
25551 avg_incubation_days,
25559 template<
typename TSeq>
25561 std::vector< double > proportions_,
25567 create_init_function_seir<TSeq>(proportions_)
25574 template<
typename TSeq>
25582 std::vector< int > h_date;
25583 std::vector< std::string > h_state;
25584 std::vector< int > h_counts;
25586 this_const->get_db().get_hist_total(
25593 std::vector< double > S(this_const->get_ndays(), 0.0);
25594 for (
size_t i = 0; i < h_date.size(); ++i)
25596 if (h_state[i] ==
"Susceptible")
25597 S[h_date[i]] += h_counts[i];
25601 double days_exposed = this_const->par(
"Avg. Incubation days");
25606 std::vector< double > gen_times(
25607 this_const->get_ndays(), 1.0 + days_exposed
25610 double p_c = this_const->par(
"Contact rate")/this_const->size();
25611 double p_i = this_const->par(
"Prob. Transmission");
25612 double p_r = this_const->par(
"Prob. Recovery");
25614 for (
size_t i = 0u; i < this_const->get_ndays(); ++i)
25616 gen_times[i] += gen_int_mean(
25774 create_init_function_sird<TSeq>(proportions_)
26137 create_init_function_seird<TSeq>(proportions_)
26329 if (which ==
static_cast<int>(m->size()))
26333 if (which ==
static_cast<int>(p->
get_id()))
26341 const auto & v = neighbor.get_virus();
26344 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
26345 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
26349 m->array_double_tmp[nviruses_tmp] =
26350 (1.0 - p->get_susceptibility_reduction(v, m)) *
26351 v->get_prob_infecting(m) *
26352 (1.0 - neighbor.get_transmission_reduction(v, m))
26355 m->array_virus_tmp[nviruses_tmp++] = &(*v);
26361 if (nviruses_tmp == 0u)
26365 int which = roulette(nviruses_tmp, m);
26370 p->set_virus(*m->array_virus_tmp[which], m);
26377 epiworld::UpdateFun<TSeq> update_infected = [](
26381 auto state = p->get_state();
26388 epiworld_fast_uint n_events = 0u;
26389 const auto & v = p->get_virus();
26392 m->array_double_tmp[n_events++] =
26393 v->get_prob_death(m) * (1.0 - p->get_death_reduction(v, m));
26396 m->array_double_tmp[n_events++] =
26397 1.0 - (1.0 - v->get_prob_recovery(m)) * (1.0 - p->get_recovery_enhancer(v, m));
26400 if (n_events == 0u)
26403 "[epi-debug] agent %i has 0 possible events!!\n",
26404 static_cast<int>(p->
get_id())
26406 throw std::logic_error(
"Zero events in exposed.");
26409 if (n_events == 0u)
26415 int which = roulette(n_events, m);
26421 if ((which % 2) == 0)
26435 throw std::logic_error(
"This function can only be applied to infected individuals. (SIR)") ;
26442 model.add_state(
"Susceptible", update_susceptible);
26443 model.add_state(
"Infected", update_infected);
26444 model.add_state(
"Recovered");
26445 model.add_state(
"Deceased");
26449 model.add_param(contact_rate,
"Contact rate");
26450 model.add_param(transmission_rate,
"Transmission rate");
26451 model.add_param(recovery_rate,
"Recovery rate");
26452 model.add_param(death_rate,
"Death rate");
26457 virus.set_state(1, 2, 3);
26458 virus.set_prob_infecting(&model(
"Transmission rate"));
26459 virus.set_prob_recovery(&model(
"Recovery rate"));
26460 virus.set_prob_death(&model(
"Death rate"));
26462 model.add_virus(virus);
26464 model.queuing_off();
26466 model.agents_empty_graph(n);
26468 model.set_name(
"Susceptible-Infected-Removed-Deceased (SIRD) (connected)");
26474 template<
typename TSeq>
26476 const std::string & vname,
26477 epiworld_fast_uint n,
26478 epiworld_double prevalence,
26479 epiworld_double contact_rate,
26480 epiworld_double transmission_rate,
26481 epiworld_double recovery_rate,
26482 epiworld_double death_rate
26710 if (which ==
static_cast<int>(ninfected))
26720 const auto & v = neighbor.get_virus();
26723 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
26724 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
26728 m->array_double_tmp[nviruses_tmp] =
26729 (1.0 - p->get_susceptibility_reduction(v, m)) *
26730 v->get_prob_infecting(m) *
26731 (1.0 - neighbor.get_transmission_reduction(v, m))
26734 m->array_virus_tmp[nviruses_tmp++] = &(*v);
26738 if (nviruses_tmp == 0u)
26742 int which = roulette(nviruses_tmp, m);
26748 *m->array_virus_tmp[which],
26757 epiworld::UpdateFun<TSeq> update_infected = [](
26761 auto state = p->get_state();
26767 auto & v = p->get_virus();
26770 if (m->runif() < 1.0/(v->get_incubation(m)))
26783 epiworld_fast_uint n_events = 0u;
26784 const auto & v = p->get_virus();
26787 m->array_double_tmp[n_events++] =
26788 v->get_prob_death(m) * (1.0 - p->get_death_reduction(v, m));
26791 m->array_double_tmp[n_events++] =
26792 1.0 - (1.0 - v->get_prob_recovery(m)) * (1.0 - p->get_recovery_enhancer(v, m));
26795 if (n_events == 0u)
26798 "[epi-debug] agent %i has 0 possible events!!\n",
26799 static_cast<int>(p->
get_id())
26801 throw std::logic_error(
"Zero events in exposed.");
26804 if (n_events == 0u)
26810 int which = roulette(n_events, m);
26816 if ((which % 2) == 0)
26830 throw std::logic_error(
"This function can only be applied to exposed or infected individuals. (SEIRD)") ;
26837 model.add_param(contact_rate,
"Contact rate");
26838 model.add_param(transmission_rate,
"Prob. Transmission");
26839 model.add_param(recovery_rate,
"Prob. Recovery");
26840 model.add_param(avg_incubation_days,
"Avg. Incubation days");
26841 model.add_param(death_rate,
"Death rate");
26844 model.add_state(
"Susceptible", update_susceptible);
26845 model.add_state(
"Exposed", update_infected);
26846 model.add_state(
"Infected", update_infected);
26847 model.add_state(
"Removed");
26848 model.add_state(
"Deceased");
26856 if (model ==
nullptr)
26857 throw std::logic_error(
26858 std::string(
"Internal error in the ModelSEIRDCONN model: ") +
26859 std::string(
"The model returns a null pointer.")
26862 model->update_infected();
26867 model.add_globalevent(update,
"Update infected individuals");
26878 virus.set_prob_infecting(&model(
"Prob. Transmission"));
26879 virus.set_prob_recovery(&model(
"Prob. Recovery"));
26880 virus.set_incubation(&model(
"Avg. Incubation days"));
26881 virus.set_prob_death(&model(
"Death rate"));
26882 model.add_virus(virus);
26884 model.queuing_off();
26887 model.agents_empty_graph(n);
26889 model.set_name(
"Susceptible-Exposed-Infected-Removed-Deceased (SEIRD) (connected)");
26895 template<
typename TSeq>
26897 const std::string & vname,
26898 epiworld_fast_uint n,
26899 epiworld_double prevalence,
26900 epiworld_double contact_rate,
26901 epiworld_double transmission_rate,
26902 epiworld_double avg_incubation_days,
26903 epiworld_double recovery_rate,
26904 epiworld_double death_rate
26915 avg_incubation_days,
26924 template<
typename TSeq>
26926 std::vector< double > proportions_,
26931 create_init_function_seird<TSeq>(proportions_)
27085 for (
const auto & c : coef_infect_cols)
27088 throw std::range_error(
"Columns specified in coef_infect_cols out of range.");
27091 for (
const auto & c : coef_recover_cols)
27094 throw std::range_error(
"Columns specified in coef_recover_cols out of range.");
27098 if (coefs_infect.size() != (coef_infect_cols.size() + 1u))
27099 throw std::logic_error(
27100 "The number of coefficients (infection) doesn't match the number of features. It must be as many features of the agents plus 1 (exposure.)"
27103 if (coefs_recover.size() != coef_recover_cols.size())
27104 throw std::logic_error(
27105 "The number of coefficients (recovery) doesn't match the number of features. It must be as many features of the agents."
27124 template<
typename TSeq>
27127 const std::string & vname,
27130 std::vector< double > coefs_infect,
27131 std::vector< double > coefs_recover,
27132 std::vector< size_t > coef_infect_cols,
27133 std::vector< size_t > coef_recover_cols,
27134 epiworld_double transmission_rate,
27135 epiworld_double recovery_rate,
27136 epiworld_double prevalence
27140 if (coef_infect_cols.size() == 0u)
27141 throw std::logic_error(
"No columns specified for coef_infect_cols.");
27143 if (coef_recover_cols.size() == 0u)
27144 throw std::logic_error(
"No columns specified for coef_recover_cols.");
27147 model.set_agents_data(
27151 model.coefs_infect = coefs_infect;
27152 model.coefs_recover = coefs_recover;
27153 model.coef_infect_cols = coef_infect_cols;
27154 model.coef_recover_cols = coef_recover_cols;
27156 epiworld::UpdateFun<TSeq> update_susceptible = [](
27165 const double coef_exposure = _m->coefs_infect[0u];
27168 size_t nviruses_tmp = 0u;
27170 double baseline = 0.0;
27171 for (
size_t k = 0u; k < _m->coef_infect_cols.size(); ++k)
27172 baseline += p->operator[](k) * _m->coefs_infect[k + 1u];
27174 for (
auto & neighbor: p->get_neighbors())
27177 if (neighbor->get_virus() ==
nullptr)
27180 auto & v = neighbor->get_virus();
27183 if (nviruses_tmp >= m->array_virus_tmp.size())
27184 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
27188 m->array_double_tmp[nviruses_tmp] =
27190 (1.0 - p->get_susceptibility_reduction(v, m)) *
27191 v->get_prob_infecting(m) *
27192 (1.0 - neighbor->get_transmission_reduction(v, m)) *
27197 m->array_double_tmp[nviruses_tmp] = 1.0/
27198 (1.0 + std::exp(-m->array_double_tmp[nviruses_tmp]));
27200 m->array_virus_tmp[nviruses_tmp++] = &(*v);
27205 if (nviruses_tmp == 0u)
27209 int which = roulette(nviruses_tmp, m);
27214 p->set_virus(*m->array_virus_tmp[which], m);
27220 epiworld::UpdateFun<TSeq> update_infected = [](
27230 #if defined(__OPENMP) || defined(_OPENMP)
27231 #pragma omp simd reduction(+:prob)
27233 for (
size_t i = 0u; i < _m->coefs_recover.size(); ++i)
27234 prob += p->operator[](i) * _m->coefs_recover[i];
27237 prob = 1.0/(1.0 + std::exp(-prob));
27239 if (prob > m->runif())
27247 model.add_state(
"Susceptible", update_susceptible);
27248 model.add_state(
"Infected", update_infected);
27249 model.add_state(
"Recovered");
27253 model.add_param(transmission_rate,
"Transmission rate");
27254 model.add_param(recovery_rate,
"Recovery rate");
27265 virus.set_prob_infecting(&model(
"Transmission rate"));
27266 virus.set_prob_recovery(&model(
"Recovery rate"));
27270 model.add_virus(virus);
27272 model.set_name(
"Susceptible-Infected-Removed (SIR) (logit)");
27278 template<
typename TSeq>
27280 const std::string & vname,
27283 std::vector< double > coefs_infect,
27284 std::vector< double > coefs_recover,
27285 std::vector< size_t > coef_infect_cols,
27286 std::vector< size_t > coef_recover_cols,
27287 epiworld_double transmission_rate,
27288 epiworld_double recovery_rate,
27289 epiworld_double prevalence
27431 (1.0 - agent.get_susceptibility_reduction(v, m)) *
27432 (1.0 - agent.get_transmission_reduction(v, m))
27435 size_t vid = v->get_id();
27438 stored[vid] =
true;
27439 innovations[vid] = &(*v);
27441 exposure[vid] += p_i;
27449 for (
size_t i = 0u; i < nviruses; ++i)
27452 if (diffmodel->normalize_exposure)
27453 exposure.at(i) /= agent.get_n_neighbors();
27455 for (
auto & j: diffmodel->data_cols)
27456 exposure.at(i) += agent(j) * diffmodel->params.at(j);
27459 double p = m->get_viruses()[i]->get_prob_infecting(m);
27460 exposure.at(i) += std::log(p) - std::log(1.0 - p);
27463 exposure.at(i) = 1.0/(1.0 + std::exp(-exposure.at(i)));
27468 int which = roulette<int>(exposure, m);
27476 *innovations.at(which),
27478 ModelDiffNet::ADOPTER
27486 model.set_agents_data(agents_data, data_ncols);
27489 model.add_state(
"Non adopters", update_non_adopters);
27490 model.add_state(
"Adopters");
27493 std::string parname = std::string(
"Prob. Adopting ") + innovation_name;
27494 model.add_param(prob_adopt, parname);
27498 innovation.set_state(1,1,1);
27500 innovation.set_prob_infecting(&model(parname));
27502 model.add_virus(innovation);
27505 std::string(
"Diffusion of Innovations - ") + innovation_name);
27511 template<
typename TSeq>
27513 const std::string & innovation_name,
27514 epiworld_double prevalence,
27515 epiworld_double prob_adopt,
27516 bool normalize_exposure,
27517 double * agents_data,
27519 std::vector< size_t > data_cols,
27520 std::vector< double > params
27529 normalize_exposure,
27571 [[assume((output) !=
nullptr)]];
27574 #define GET_MODEL(model, output) \
27575 auto * output = dynamic_cast< ModelSEIRMixing<TSeq> * >( (model) ); \
27576 assert((output) != nullptr);
27583 template<
typename TSeq = EPI_DEFAULT_TSEQ>
27589 std::vector< size_t > infected;
27592 std::vector< size_t > n_infected_per_group;
27595 std::vector< size_t > entity_indices;
27597 void update_infected_list();
27598 std::vector< size_t > sampled_agents;
27599 size_t sample_agents(
27601 std::vector< size_t > & sampled_agents
27603 std::vector< double > adjusted_contact_rate;
27604 std::vector< double > contact_matrix;
27607 std::vector< int > sampled_sizes;
27612 static const int SUSCEPTIBLE = 0;
27613 static const int EXPOSED = 1;
27614 static const int INFECTED = 2;
27615 static const int RECOVERED = 3;
27635 const std::string & vname,
27636 epiworld_fast_uint n,
27637 epiworld_double prevalence,
27638 epiworld_double contact_rate,
27639 epiworld_double transmission_rate,
27640 epiworld_double avg_incubation_days,
27641 epiworld_double recovery_rate,
27642 std::vector< double > contact_matrix
27658 const std::string & vname,
27659 epiworld_fast_uint n,
27660 epiworld_double prevalence,
27661 epiworld_double contact_rate,
27662 epiworld_double transmission_rate,
27663 epiworld_double avg_incubation_days,
27664 epiworld_double recovery_rate,
27665 std::vector< double > contact_matrix
27669 epiworld_fast_uint ndays,
27683 std::vector< double > proportions_,
27684 std::vector< int > queue_ = {}
27687 void set_contact_matrix(std::vector< double > cmat)
27689 contact_matrix = cmat;
27695 template<
typename TSeq>
27701 std::fill(n_infected_per_group.begin(), n_infected_per_group.end(), 0u);
27703 for (
auto & a : agents)
27708 if (a.get_n_entities() > 0u)
27710 const auto & entity = a.get_entity(0u);
27713 entity_indices[entity.get_id()] +
27715 n_infected_per_group[entity.get_id()]++
27727 template<
typename TSeq>
27730 std::vector< size_t > & sampled_agents
27734 size_t agent_group_id = agent->get_entity(0u).get_id();
27735 size_t ngroups = this->entities.size();
27738 for (
size_t g = 0; g < ngroups; ++g)
27741 size_t group_size = n_infected_per_group[g];
27746 adjusted_contact_rate[g] * contact_matrix[
27747 MM(agent_group_id, g, ngroups)
27755 for (
int s = 0; s < nsamples; ++s)
27762 if (which >=
static_cast<int>(group_size))
27763 which =
static_cast<int>(group_size) - 1;
27766 auto & a = this->population.at(infected.at(entity_indices[g] + which));
27768 auto & a = this->get_agent(infected[entity_indices[g] + which]);
27773 throw std::logic_error(
27774 "The agent is not infected, but it should be."
27779 if (a.get_id() == agent->
get_id())
27782 sampled_agents[samp_id++] = a.get_id();
27792 template<
typename TSeq>
27794 epiworld_fast_uint ndays,
27805 template<
typename TSeq>
27812 size_t nentities = this->entities.size();
27813 if (this->contact_matrix.size() != nentities*nentities)
27814 throw std::length_error(
27815 std::string(
"The contact matrix must be a square matrix of size ") +
27816 std::string(
"nentities x nentities. ") +
27817 std::to_string(this->contact_matrix.size()) +
27818 std::string(
" != ") + std::to_string(nentities*nentities) +
27822 for (
size_t i = 0u; i < this->entities.size(); ++i)
27825 for (
size_t j = 0u; j < this->entities.size(); ++j)
27827 if (this->contact_matrix[MM(i, j, nentities)] < 0.0)
27828 throw std::range_error(
27829 std::string(
"The contact matrix must be non-negative. ") +
27830 std::to_string(this->contact_matrix[MM(i, j, nentities)]) +
27831 std::string(
" < 0.")
27833 sum += this->contact_matrix[MM(i, j, nentities)];
27835 if (sum < 0.999 || sum > 1.001)
27836 throw std::range_error(
27837 std::string(
"The contact matrix must have rows that add to one. ") +
27838 std::to_string(sum) +
27839 std::string(
" != 1.")
27847 n_infected_per_group.resize(this->entities.size(), 0u);
27848 std::fill(n_infected_per_group.begin(), n_infected_per_group.end(), 0u);
27852 std::fill(infected.begin(), infected.end(), 0u);
27855 entity_indices.resize(this->entities.size(), 0u);
27856 std::fill(entity_indices.begin(), entity_indices.end(), 0u);
27857 for (
size_t i = 1u; i < this->entities.size(); ++i)
27860 entity_indices[i] +=
27861 this->entities[i - 1].size() +
27862 entity_indices[i - 1]
27868 adjusted_contact_rate.clear();
27869 adjusted_contact_rate.resize(this->entities.size(), 0.0);
27871 for (
size_t i = 0u; i < this->entities.size(); ++i)
27874 adjusted_contact_rate[i] =
27876 static_cast< epiworld_double
> (this->get_entity(i).size());
27880 if (adjusted_contact_rate[i] > 1.0)
27881 adjusted_contact_rate[i] = 1.0;
27885 this->update_infected_list();
27891 template<
typename TSeq>
27899 #if __cplusplus >= 202302L
27901 [[assume(ptr !=
nullptr)]];
27904 assert(ptr !=
nullptr);
27922 template<
typename TSeq>
27925 const std::string & vname,
27926 epiworld_fast_uint n,
27927 epiworld_double prevalence,
27928 epiworld_double contact_rate,
27929 epiworld_double transmission_rate,
27930 epiworld_double avg_incubation_days,
27931 epiworld_double recovery_rate,
27932 std::vector< double > contact_matrix
27937 this->contact_matrix = contact_matrix;
27939 epiworld::UpdateFun<TSeq> update_susceptible = [](
27944 if (p->get_n_entities() == 0)
27949 GET_MODEL(m, m_down);
27951 size_t ndraws = m_down->sample_agents(p, m_down->sampled_agents);
27954 m_down->sampled_sizes.push_back(
static_cast<int>(ndraws));
27961 int nviruses_tmp = 0;
27962 for (
size_t n = 0u; n < ndraws; ++n)
27965 auto & neighbor = m->get_agent(m_down->sampled_agents[n]);
27967 auto & v = neighbor.get_virus();
27970 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
27971 throw std::logic_error(
27972 "Trying to add an extra element to a temporal array outside of the range."
27977 m->array_double_tmp[nviruses_tmp] =
27978 (1.0 - p->get_susceptibility_reduction(v, m)) *
27979 v->get_prob_infecting(m) *
27980 (1.0 - neighbor.get_transmission_reduction(v, m))
27983 m->array_virus_tmp[nviruses_tmp++] = &(*v);
27988 int which = roulette(nviruses_tmp, m);
27994 *m->array_virus_tmp[which],
28003 epiworld::UpdateFun<TSeq> update_exposed_and_infected = [](
28007 auto state = p->get_state();
28013 auto & v = p->get_virus();
28016 if (m->runif() < 1.0/(v->get_incubation(m)))
28030 epiworld_fast_uint n_events = 0u;
28031 const auto & v = p->get_virus();
28034 m->array_double_tmp[n_events++] =
28035 1.0 - (1.0 - v->get_prob_recovery(m)) * (1.0 - p->get_recovery_enhancer(v, m));
28038 if (n_events == 0u)
28041 "[epi-debug] agent %i has 0 possible events!!\n",
28042 static_cast<int>(p->
get_id())
28044 throw std::logic_error(
"Zero events in exposed.");
28047 if (n_events == 0u)
28053 int which = roulette(n_events, m);
28064 throw std::logic_error(
"This function can only be applied to exposed or infected individuals. (SEIR)") ;
28071 model.add_param(contact_rate,
"Contact rate");
28072 model.add_param(transmission_rate,
"Prob. Transmission");
28073 model.add_param(recovery_rate,
"Prob. Recovery");
28074 model.add_param(avg_incubation_days,
"Avg. Incubation days");
28077 model.add_state(
"Susceptible", update_susceptible);
28078 model.add_state(
"Exposed", update_exposed_and_infected);
28079 model.add_state(
"Infected", update_exposed_and_infected);
28080 model.add_state(
"Recovered");
28086 GET_MODEL(m, m_down);
28088 m_down->update_infected_list();
28094 model.add_globalevent(update,
"Update infected individuals");
28105 virus.set_prob_infecting(&model(
"Prob. Transmission"));
28106 virus.set_prob_recovery(&model(
"Prob. Recovery"));
28107 virus.set_incubation(&model(
"Avg. Incubation days"));
28109 model.add_virus(virus);
28111 model.queuing_off();
28114 model.agents_empty_graph(n);
28116 model.set_name(
"Susceptible-Exposed-Infected-Removed (SEIR) with Mixing");
28122 template<
typename TSeq>
28124 const std::string & vname,
28125 epiworld_fast_uint n,
28126 epiworld_double prevalence,
28127 epiworld_double contact_rate,
28128 epiworld_double transmission_rate,
28129 epiworld_double avg_incubation_days,
28130 epiworld_double recovery_rate,
28131 std::vector< double > contact_matrix
28135 this->contact_matrix = contact_matrix;
28144 avg_incubation_days,
28153 template<
typename TSeq>
28155 std::vector< double > proportions_,
28161 create_init_function_seir<TSeq>(proportions_)
28197 [[assume((output) !=
nullptr)]];
28200 #define GET_MODEL(model, output) \
28201 auto * output = dynamic_cast< ModelSIRMixing<TSeq> * >( (model) ); \
28202 assert((output) != nullptr);
28209 template<
typename TSeq = EPI_DEFAULT_TSEQ>
28215 std::vector< size_t > infected;
28219 std::vector< size_t > n_infected_per_group;
28222 std::vector< size_t > entity_indices;
28224 void update_infected_list();
28225 std::vector< size_t > sampled_agents;
28226 size_t sample_agents(
28228 std::vector< size_t > & sampled_agents
28231 std::vector< double > adjusted_contact_rate;
28232 std::vector< double > contact_matrix;
28234 size_t index(
size_t i,
size_t j,
size_t n) {
28240 static const int SUSCEPTIBLE = 0;
28241 static const int INFECTED = 1;
28242 static const int RECOVERED = 2;
28260 const std::string & vname,
28261 epiworld_fast_uint n,
28262 epiworld_double prevalence,
28263 epiworld_double contact_rate,
28264 epiworld_double transmission_rate,
28265 epiworld_double recovery_rate,
28266 std::vector< double > contact_matrix
28281 const std::string & vname,
28282 epiworld_fast_uint n,
28283 epiworld_double prevalence,
28284 epiworld_double contact_rate,
28285 epiworld_double transmission_rate,
28286 epiworld_double recovery_rate,
28287 std::vector< double > contact_matrix
28291 epiworld_fast_uint ndays,
28305 std::vector< double > proportions_,
28306 std::vector< int > queue_ = {}
28309 size_t get_n_infected(
size_t group)
const
28311 return n_infected_per_group[group];
28314 void set_contact_matrix(std::vector< double > cmat)
28316 contact_matrix = cmat;
28322 template<
typename TSeq>
28327 std::fill(n_infected_per_group.begin(), n_infected_per_group.end(), 0u);
28330 for (
auto & a : agents)
28334 if (a.get_n_entities() > 0u)
28336 const auto & entity = a.get_entity(0u);
28339 entity_indices[entity.get_id()] +
28341 n_infected_per_group[entity.get_id()]++
28353 template<
typename TSeq>
28356 std::vector< size_t > & sampled_agents
28360 size_t agent_group_id = agent->get_entity(0u).get_id();
28361 size_t ngroups = this->entities.size();
28364 for (
size_t g = 0; g < ngroups; ++g)
28367 size_t group_size = n_infected_per_group[g];
28372 adjusted_contact_rate[g] * contact_matrix[
28373 index(agent_group_id, g, ngroups)
28381 for (
int s = 0; s < nsamples; ++s)
28388 if (which >=
static_cast<int>(group_size))
28389 which =
static_cast<int>(group_size) - 1;
28392 auto & a = this->population.at(infected.at(entity_indices[g] + which));
28394 auto & a = this->get_agent(infected[entity_indices[g] + which]);
28399 throw std::logic_error(
28400 "The agent is not infected, but it should be."
28405 if (a.get_id() == agent->
get_id())
28408 sampled_agents[samp_id++] = a.get_id();
28418 template<
typename TSeq>
28420 epiworld_fast_uint ndays,
28430 template<
typename TSeq>
28437 size_t nentities = this->entities.size();
28438 if (this->contact_matrix.size() != nentities*nentities)
28439 throw std::length_error(
28440 std::string(
"The contact matrix must be a square matrix of size ") +
28441 std::string(
"nentities x nentities. ") +
28442 std::to_string(this->contact_matrix.size()) +
28443 std::string(
" != ") + std::to_string(nentities*nentities) +
28447 for (
size_t i = 0u; i < this->entities.size(); ++i)
28450 for (
size_t j = 0u; j < this->entities.size(); ++j)
28452 if (this->contact_matrix[index(i, j, nentities)] < 0.0)
28453 throw std::range_error(
28454 std::string(
"The contact matrix must be non-negative. ") +
28455 std::to_string(this->contact_matrix[index(i, j, nentities)]) +
28456 std::string(
" < 0.")
28458 sum += this->contact_matrix[index(i, j, nentities)];
28460 if (sum < 0.999 || sum > 1.001)
28461 throw std::range_error(
28462 std::string(
"The contact matrix must have rows that add to one. ") +
28463 std::to_string(sum) +
28464 std::string(
" != 1.")
28472 n_infected_per_group.resize(this->entities.size(), 0u);
28473 std::fill(n_infected_per_group.begin(), n_infected_per_group.end(), 0u);
28477 std::fill(infected.begin(), infected.end(), 0u);
28480 entity_indices.resize(this->entities.size(), 0u);
28481 std::fill(entity_indices.begin(), entity_indices.end(), 0u);
28482 for (
size_t i = 1u; i < this->entities.size(); ++i)
28484 entity_indices[i] +=
28485 this->entities[i - 1].size() +
28486 entity_indices[i - 1]
28491 adjusted_contact_rate.clear();
28492 adjusted_contact_rate.resize(this->entities.size(), 0.0);
28494 for (
size_t i = 0u; i < this->entities.size(); ++i)
28496 adjusted_contact_rate[i] =
28498 static_cast< epiworld_double
> (this->get_entity(i).size());
28501 if (adjusted_contact_rate[i] > 1.0)
28502 adjusted_contact_rate[i] = 1.0;
28505 this->update_infected_list();
28510 template<
typename TSeq>
28533 template<
typename TSeq>
28536 const std::string & vname,
28537 epiworld_fast_uint n,
28538 epiworld_double prevalence,
28539 epiworld_double contact_rate,
28540 epiworld_double transmission_rate,
28541 epiworld_double recovery_rate,
28542 std::vector< double > contact_matrix
28547 this->contact_matrix = contact_matrix;
28549 epiworld::UpdateFun<TSeq> update_susceptible = [](
28554 if (p->get_n_entities() == 0)
28559 GET_MODEL(m, m_down);
28561 size_t ndraws = m_down->sample_agents(p, m_down->sampled_agents);
28568 int nviruses_tmp = 0;
28569 for (
size_t n = 0u; n < ndraws; ++n)
28572 auto & neighbor = m->get_agent(m_down->sampled_agents[n]);
28574 auto & v = neighbor.get_virus();
28577 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
28578 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
28582 m->array_double_tmp[nviruses_tmp] =
28583 (1.0 - p->get_susceptibility_reduction(v, m)) *
28584 v->get_prob_infecting(m) *
28585 (1.0 - neighbor.get_transmission_reduction(v, m))
28588 m->array_virus_tmp[nviruses_tmp++] = &(*v);
28593 int which = roulette(nviruses_tmp, m);
28599 *m->array_virus_tmp[which],
28608 epiworld::UpdateFun<TSeq> update_infected = [](
28612 auto state = p->get_state();
28619 epiworld_fast_uint n_events = 0u;
28620 const auto & v = p->get_virus();
28623 m->array_double_tmp[n_events++] =
28624 1.0 - (1.0 - v->get_prob_recovery(m)) * (1.0 - p->get_recovery_enhancer(v, m));
28627 if (n_events == 0u)
28630 "[epi-debug] agent %i has 0 possible events!!\n",
28631 static_cast<int>(p->
get_id())
28633 throw std::logic_error(
"Zero events in infected.");
28636 if (n_events == 0u)
28642 int which = roulette(n_events, m);
28653 throw std::logic_error(
"This function can only be applied to infected individuals. (SIR)") ;
28660 model.add_param(contact_rate,
"Contact rate");
28661 model.add_param(transmission_rate,
"Prob. Transmission");
28662 model.add_param(recovery_rate,
"Prob. Recovery");
28665 model.add_state(
"Susceptible", update_susceptible);
28666 model.add_state(
"Infected", update_infected);
28667 model.add_state(
"Recovered");
28673 GET_MODEL(m, m_down);
28675 m_down->update_infected_list();
28681 model.add_globalevent(update,
"Update infected individuals");
28692 virus.set_prob_infecting(&model(
"Prob. Transmission"));
28693 virus.set_prob_recovery(&model(
"Prob. Recovery"));
28695 model.add_virus(virus);
28697 model.queuing_off();
28700 model.agents_empty_graph(n);
28702 model.set_name(
"Susceptible-Infected-Removed (SIR) with Mixing");
28708 template<
typename TSeq>
28710 const std::string & vname,
28711 epiworld_fast_uint n,
28712 epiworld_double prevalence,
28713 epiworld_double contact_rate,
28714 epiworld_double transmission_rate,
28715 epiworld_double recovery_rate,
28716 std::vector< double > contact_matrix
28720 this->contact_matrix = contact_matrix;
28737 template<
typename TSeq>
28739 std::vector< double > proportions_,
28745 create_init_function_sir<TSeq>(proportions_)
29140 if (which ==
static_cast<int>(n_infectious))
29146 if (neighbor.
get_id() == p->get_id())
29153 if (neighbor.get_virus() ==
nullptr)
29154 throw std::logic_error(
"The neighbor has no virus.");
29157 if (neighbor.get_state() != model->PRODROMAL)
29158 throw std::logic_error(
29159 "The neighbor is not in the prodromal state. The state is: " +
29160 std::to_string(neighbor.get_state())
29163 auto & v = neighbor.get_virus();
29166 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
29167 throw std::logic_error(
"Trying to add an extra element to a temporal array outside of the range.");
29171 m->array_double_tmp[nviruses_tmp] =
29172 (1.0 - p->get_susceptibility_reduction(v, m)) *
29173 v->get_prob_infecting(m) *
29174 (1.0 - neighbor.get_transmission_reduction(v, m))
29177 m->array_virus_tmp[nviruses_tmp++] = &(*v);
29182 if (nviruses_tmp == 0u)
29186 int which = roulette(nviruses_tmp, m);
29191 p->set_virus(*m->array_virus_tmp[which], m);
29197 LOCAL_UPDATE_FUN(m_update_exposed) {
29199 if (m->runif() < (1.0/p->get_virus()->get_incubation(m)))
29206 LOCAL_UPDATE_FUN(m_update_prodromal) {
29208 if (m->runif() < (1.0/m->par(
"Prodromal period")))
29211 GET_MODEL(m, model);
29212 model->day_rash_onset[p->get_id()] = m->
today();
29221 LOCAL_UPDATE_FUN(m_update_rash) {
29224 GET_MODEL(m, model);
29227 if (
static_cast<int>(model->day_flagged.size()) <= p->get_id())
29228 throw std::logic_error(
29229 "The agent is not in the list of quarantined or isolated agents: " +
29230 std::to_string(p->get_id()) +
29232 std::to_string(model->day_flagged.size()) +
29233 ". The model has " + std::to_string(model->size()) +
" agents."
29239 bool detected =
false;
29241 (m->par(
"Isolation period") >= 0) &&
29242 (m->runif() < 1.0/m->par(
"Days undetected"))
29245 model->system_quarantine_triggered =
true;
29252 m->array_double_tmp[0] = 1.0/m->par(
"Rash period");
29253 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
29264 ModelMeaslesSchool::ISOLATED_RECOVERED:
29265 ModelMeaslesSchool::RECOVERED
29268 else if (which == 1)
29275 ModelMeaslesSchool::DETECTED_HOSPITALIZED :
29276 ModelMeaslesSchool::HOSPITALIZED
29279 else if (which != 0)
29281 throw std::logic_error(
"The roulette returned an unexpected value.");
29283 else if ((which == 0u) && detected)
29287 p->change_state(m, ModelMeaslesSchool::ISOLATED);
29292 LOCAL_UPDATE_FUN(m_update_isolated) {
29294 GET_MODEL(m, model);
29298 int days_since = m->
today() - model->day_rash_onset[p->get_id()];
29301 (m->par(
"Isolation period") <= days_since) ?
29306 m->array_double_tmp[0] = 1.0/m->par(
"Rash period");
29307 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
29319 ModelMeaslesSchool::RECOVERED
29324 m, ModelMeaslesSchool::ISOLATED_RECOVERED
29329 else if (which == 1u)
29335 ModelMeaslesSchool::HOSPITALIZED :
29336 ModelMeaslesSchool::DETECTED_HOSPITALIZED
29341 else if ((which == 0u) && unisolate)
29343 p->change_state(m, ModelMeaslesSchool::RASH);
29348 LOCAL_UPDATE_FUN(m_update_isolated_recovered) {
29350 GET_MODEL(m, model);
29354 int days_since = m->
today() - model->day_rash_onset[p->get_id()];
29357 (m->par(
"Isolation period") <= days_since) ?
29361 p->change_state(m, ModelMeaslesSchool::RECOVERED);
29365 LOCAL_UPDATE_FUN(m_update_q_exposed) {
29368 GET_MODEL(m, model);
29370 m->
today() - model->day_flagged[p->get_id()];
29372 bool unquarantine =
29373 (m->par(
"Quarantine period") <= days_since) ?
29377 if (m->runif() < (1.0/p->get_virus()->get_incubation(m)))
29385 ModelMeaslesSchool::PRODROMAL
29390 ModelMeaslesSchool::QUARANTINED_PRODROMAL
29394 else if (unquarantine)
29398 ModelMeaslesSchool::EXPOSED
29404 LOCAL_UPDATE_FUN(m_update_q_susceptible) {
29406 GET_MODEL(m, model);
29408 m->
today() - model->day_flagged[p->get_id()];
29410 if (days_since >= m->par(
"Quarantine period"))
29411 p->change_state(m, ModelMeaslesSchool::SUSCEPTIBLE);
29415 LOCAL_UPDATE_FUN(m_update_q_prodromal) {
29417 GET_MODEL(m, model);
29421 int days_since = m->
today() - model->day_flagged[p->get_id()];
29423 bool unquarantine =
29424 (m->par(
"Quarantine period") <= days_since) ?
29428 if (m->runif() < (1.0/m->par(
"Prodromal period")))
29430 model->day_rash_onset[p->get_id()] = m->
today();
29431 p->change_state(m, ModelMeaslesSchool::ISOLATED);
29437 p->change_state(m, ModelMeaslesSchool::PRODROMAL);
29443 LOCAL_UPDATE_FUN(m_update_q_recovered) {
29445 GET_MODEL(m, model);
29446 int days_since = m->
today() - model->day_flagged[p->get_id()];
29448 if (days_since >= m->par(
"Quarantine period"))
29449 p->change_state(m, ModelMeaslesSchool::RECOVERED);
29453 LOCAL_UPDATE_FUN(m_update_hospitalized) {
29456 if (m->runif() < 1.0/m->par(
"Hospitalization period"))
29457 p->rm_virus(m, ModelMeaslesSchool::RECOVERED);
29464 template<
typename TSeq>
29467 epiworld_fast_uint n,
29468 epiworld_fast_uint n_exposed,
29470 epiworld_double contact_rate,
29471 epiworld_double transmission_rate,
29472 epiworld_double vax_efficacy,
29473 epiworld_double vax_reduction_recovery_rate,
29474 epiworld_double incubation_period,
29475 epiworld_double prodromal_period,
29476 epiworld_double rash_period,
29477 epiworld_double days_undetected,
29478 epiworld_double hospitalization_rate,
29479 epiworld_double hospitalization_period,
29481 epiworld_double prop_vaccinated,
29482 epiworld_fast_int quarantine_period,
29483 epiworld_double quarantine_willingness,
29484 epiworld_fast_int isolation_period
29487 model.add_state(
"Susceptible", this->m_update_susceptible);
29488 model.add_state(
"Exposed", this->m_update_exposed);
29489 model.add_state(
"Prodromal", this->m_update_prodromal);
29490 model.add_state(
"Rash", this->m_update_rash);
29491 model.add_state(
"Isolated", this->m_update_isolated);
29493 "Isolated Recovered", this->m_update_isolated_recovered
29495 model.add_state(
"Detected Hospitalized", this->m_update_hospitalized);
29497 "Quarantined Exposed", this->m_update_q_exposed
29501 "Quarantined Susceptible", this->m_update_q_susceptible
29505 "Quarantined Prodromal", this->m_update_q_prodromal
29509 "Quarantined Recovered", this->m_update_q_recovered
29512 model.add_state(
"Hospitalized", this->m_update_hospitalized);
29514 model.add_state(
"Recovered");
29517 model.add_param(contact_rate,
"Contact rate");
29518 model.add_param(transmission_rate,
"Transmission rate");
29519 model.add_param(incubation_period,
"Incubation period");
29520 model.add_param(prodromal_period,
"Prodromal period");
29521 model.add_param(rash_period,
"Rash period");
29522 model.add_param(days_undetected,
"Days undetected");
29523 model.add_param(quarantine_period,
"Quarantine period");
29525 quarantine_willingness,
"Quarantine willingness"
29527 model.add_param(isolation_period,
"Isolation period");
29528 model.add_param(hospitalization_rate,
"Hospitalization rate");
29529 model.add_param(hospitalization_period,
"Hospitalization period");
29530 model.add_param(prop_vaccinated,
"Vaccination rate");
29531 model.add_param(vax_efficacy,
"Vax efficacy");
29532 model.add_param(vax_reduction_recovery_rate,
"(IGNORED) Vax improved recovery");
29536 measles.set_state(EXPOSED, RECOVERED);
29537 measles.set_prob_infecting(&model(
"Transmission rate"));
29538 measles.set_prob_recovery(&model(
"Rash period"));
29539 measles.set_incubation(&model(
"Incubation period"));
29540 measles.set_distribution(
29541 distribute_virus_randomly(n_exposed,
false)
29544 model.add_virus(measles);
29547 Tool<> vaccine(
"Vaccine");
29548 vaccine.set_susceptibility_reduction(&model(
"Vax efficacy"));
29549 vaccine.set_recovery_enhancer(&model(
"(IGNORED) Vax improved recovery"));
29550 vaccine.set_distribution(
29551 distribute_tool_randomly(prop_vaccinated,
true)
29554 model.add_tool(vaccine);
29557 model.add_globalevent(this->m_update_model,
"Update model");
29558 model.queuing_off();
29561 model.agents_empty_graph(n);
29567 template<
typename TSeq>
29569 epiworld_fast_uint n,
29570 epiworld_fast_uint n_exposed,
29572 epiworld_double contact_rate,
29573 epiworld_double transmission_rate,
29574 epiworld_double vax_efficacy,
29575 epiworld_double vax_reduction_recovery_rate,
29576 epiworld_double incubation_period,
29577 epiworld_double prodromal_period,
29578 epiworld_double rash_period,
29579 epiworld_double days_undetected,
29580 epiworld_double hospitalization_rate,
29581 epiworld_double hospitalization_period,
29583 epiworld_double prop_vaccinated,
29584 epiworld_fast_int quarantine_period,
29585 epiworld_double quarantine_willingness,
29586 epiworld_fast_int isolation_period
29597 vax_reduction_recovery_rate,
29602 hospitalization_rate,
29603 hospitalization_period,
29606 quarantine_willingness,
29614 #undef SAMPLE_FROM_PROBS
29615 #undef LOCAL_UPDATE_FUN
29650 [[assume((output) !=
nullptr)]];
29653 #define GET_MODEL(model, output) \
29654 auto * output = dynamic_cast< ModelSEIRMixingQuarantine<TSeq> * >( (model) ); \
29655 assert((output) != nullptr);
29658 #define SAMPLE_FROM_PROBS(n, ans) \
29660 epiworld_double p_total = m->runif(); \
29661 for (ans = 0u; ans < n; ++ans) \
29663 if (p_total < m->array_double_tmp[ans]) \
29665 m->array_double_tmp[ans + 1] += m->array_double_tmp[ans]; \
29699 template<
typename TSeq = EPI_DEFAULT_TSEQ>
29705 std::vector< size_t > infected;
29708 std::vector< size_t > n_infected_per_group;
29711 std::vector< size_t > entity_indices;
29713 void m_update_infected_list();
29714 std::vector< size_t > sampled_agents;
29715 size_t sample_agents(
29717 std::vector< size_t > & sampled_agents
29719 std::vector< double > adjusted_contact_rate;
29720 std::vector< double > contact_matrix;
29723 std::vector< int > sampled_sizes;
29737 std::vector< bool > quarantine_willingness;
29738 std::vector< bool > isolation_willingness;
29739 std::vector< size_t > agent_quarantine_triggered;
29740 std::vector< int > day_flagged;
29741 std::vector< int > day_onset;
29742 std::vector< int > day_exposed;
29744 void m_quarantine_process();
29748 std::vector< size_t > tracking_matrix;
29749 std::vector< size_t > tracking_matrix_size;
29751 void m_add_tracking(
size_t infected_id,
size_t agent_id);
29755 static const int SUSCEPTIBLE = 0;
29756 static const int EXPOSED = 1;
29757 static const int INFECTED = 2;
29758 static const int ISOLATED = 3;
29759 static const int DETECTED_HOSPITALIZED = 4;
29760 static const int QUARANTINED_SUSCEPTIBLE = 5;
29761 static const int QUARANTINED_EXPOSED = 6;
29762 static const int ISOLATED_RECOVERED = 7;
29763 static const int HOSPITALIZED = 8;
29764 static const int RECOVERED = 9;
29766 static const size_t QUARANTINE_PROCESS_INACTIVE = 0u;
29767 static const size_t QUARANTINE_PROCESS_ACTIVE = 1u;
29768 static const size_t QUARANTINE_PROCESS_DONE = 2u;
29797 const std::string & vname,
29798 epiworld_fast_uint n,
29799 epiworld_double prevalence,
29800 epiworld_double contact_rate,
29801 epiworld_double transmission_rate,
29802 epiworld_double avg_incubation_days,
29803 epiworld_double recovery_rate,
29804 std::vector< double > contact_matrix,
29805 epiworld_double hospitalization_rate,
29806 epiworld_double hospitalization_period,
29808 epiworld_double days_undetected,
29809 epiworld_fast_int quarantine_period,
29810 epiworld_double quarantine_willingness,
29811 epiworld_double isolation_willingness,
29812 epiworld_fast_int isolation_period,
29813 epiworld_double contact_tracing_success_rate = 1.0,
29814 epiworld_fast_uint contact_tracing_days_prior = 4u
29839 const std::string & vname,
29840 epiworld_fast_uint n,
29841 epiworld_double prevalence,
29842 epiworld_double contact_rate,
29843 epiworld_double transmission_rate,
29844 epiworld_double avg_incubation_days,
29845 epiworld_double recovery_rate,
29846 std::vector< double > contact_matrix,
29847 epiworld_double hospitalization_rate,
29848 epiworld_double hospitalization_period,
29850 epiworld_double days_undetected,
29851 epiworld_fast_int quarantine_period,
29852 epiworld_double quarantine_willingness,
29853 epiworld_double isolation_willingness,
29854 epiworld_fast_int isolation_period,
29855 epiworld_double contact_tracing_success_rate = 1.0,
29856 epiworld_fast_uint contact_tracing_days_prior = 4u
29866 epiworld_fast_uint ndays,
29889 std::vector< double > proportions_,
29890 std::vector< int > queue_ = {}
29899 contact_matrix = cmat;
29909 return contact_matrix;
29918 return agent_quarantine_triggered;
29927 return quarantine_willingness;
29936 return isolation_willingness;
29941 template<
typename TSeq>
29944 GET_MODEL(m, model);
29945 model->m_quarantine_process();
29946 model->events_run();
29947 model->m_update_infected_list();
29951 template<
typename TSeq>
29953 size_t infected_id,
29959 if (agent_quarantine_triggered[infected_id] >= QUARANTINE_PROCESS_DONE)
29965 if (days_since_onset >
29972 size_t loc = tracking_matrix_size[infected_id] % EPI_MAX_TRACKING;
29976 tracking_matrix_size[infected_id]++;
29982 template<
typename TSeq>
29988 std::fill(n_infected_per_group.begin(), n_infected_per_group.end(), 0u);
29990 for (
auto & a : agents)
29995 if (a.get_n_entities() > 0u)
29997 const auto & entity = a.get_entity(0u);
30000 entity_indices[entity.get_id()] +
30002 n_infected_per_group[entity.get_id()]++
30014 template<
typename TSeq>
30017 std::vector< size_t > & sampled_agents
30021 size_t agent_group_id = agent->get_entity(0u).get_id();
30022 size_t ngroups = this->entities.size();
30025 for (
size_t g = 0; g < ngroups; ++g)
30028 size_t group_size = n_infected_per_group[g];
30030 if (group_size == 0u)
30036 adjusted_contact_rate[g] * contact_matrix[
30037 MM(agent_group_id, g, ngroups)
30045 for (
int s = 0; s < nsamples; ++s)
30052 if (which >=
static_cast<int>(group_size))
30053 which =
static_cast<int>(group_size) - 1;
30056 auto & a = this->population.at(infected.at(entity_indices[g] + which));
30058 auto & a = this->get_agent(infected[entity_indices[g] + which]);
30063 throw std::logic_error(
30064 "The agent is not infected, but it should be."
30069 if (a.get_id() == agent->
get_id())
30072 sampled_agents[samp_id++] = a.get_id();
30082 template<
typename TSeq>
30084 epiworld_fast_uint ndays,
30095 template<
typename TSeq>
30102 size_t nentities = this->entities.size();
30103 if (this->contact_matrix.size() != nentities*nentities)
30104 throw std::length_error(
30105 std::string(
"The contact matrix must be a square matrix of size ") +
30106 std::string(
"nentities x nentities. ") +
30107 std::to_string(this->contact_matrix.size()) +
30108 std::string(
" != ") + std::to_string(nentities*nentities) +
30112 for (
size_t i = 0u; i < this->entities.size(); ++i)
30115 for (
size_t j = 0u; j < this->entities.size(); ++j)
30117 if (this->contact_matrix[MM(i, j, nentities)] < 0.0)
30118 throw std::range_error(
30119 std::string(
"The contact matrix must be non-negative. ") +
30120 std::to_string(this->contact_matrix[MM(i, j, nentities)]) +
30121 std::string(
" < 0.")
30123 sum += this->contact_matrix[MM(i, j, nentities)];
30125 if (sum < 0.999 || sum > 1.001)
30126 throw std::range_error(
30127 std::string(
"The contact matrix must have rows that add to one. ") +
30128 std::to_string(sum) +
30129 std::string(
" != 1.")
30137 n_infected_per_group.resize(this->entities.size(), 0u);
30138 std::fill(n_infected_per_group.begin(), n_infected_per_group.end(), 0u);
30142 std::fill(infected.begin(), infected.end(), 0u);
30145 entity_indices.resize(this->entities.size(), 0u);
30146 std::fill(entity_indices.begin(), entity_indices.end(), 0u);
30147 for (
size_t i = 1u; i < this->entities.size(); ++i)
30150 entity_indices[i] +=
30151 this->entities[i - 1].size() +
30152 entity_indices[i - 1]
30158 adjusted_contact_rate.clear();
30159 adjusted_contact_rate.resize(this->entities.size(), 0.0);
30161 for (
size_t i = 0u; i < this->entities.size(); ++i)
30164 adjusted_contact_rate[i] =
30166 static_cast< epiworld_double
> (this->get_entity(i).size());
30170 if (adjusted_contact_rate[i] > 1.0)
30171 adjusted_contact_rate[i] = 1.0;
30175 this->m_update_infected_list();
30178 quarantine_willingness.resize(this->size(),
false);
30179 isolation_willingness.resize(this->size(),
false);
30180 for (
size_t idx = 0; idx < quarantine_willingness.size(); ++idx)
30182 quarantine_willingness[idx] =
30184 isolation_willingness[idx] =
30188 agent_quarantine_triggered.resize(this->size(), 0u);
30190 agent_quarantine_triggered.begin(),
30191 agent_quarantine_triggered.end(),
30195 day_flagged.resize(this->size(), 0);
30197 day_flagged.begin(),
30202 day_onset.resize(this->size(), 0);
30209 day_exposed.resize(this->size(), 0);
30211 day_exposed.begin(),
30218 std::fill(tracking_matrix.begin(), tracking_matrix.end(), 0u);
30221 std::fill(tracking_matrix_size.begin(), tracking_matrix_size.end(), 0u);
30227 template<
typename TSeq>
30235 #if __cplusplus >= 202302L
30237 [[assume(ptr !=
nullptr)]];
30240 assert(ptr !=
nullptr);
30247 template<
typename TSeq>
30252 if (p->get_n_entities() == 0)
30257 GET_MODEL(m, m_down);
30259 size_t ndraws = m_down->sample_agents(p, m_down->sampled_agents);
30262 m_down->sampled_sizes.push_back(
static_cast<int>(ndraws));
30269 int nviruses_tmp = 0;
30270 for (
size_t n = 0u; n < ndraws; ++n)
30273 auto & neighbor = m->get_agent(m_down->sampled_agents[n]);
30275 auto & v = neighbor.get_virus();
30278 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
30279 throw std::logic_error(
30280 "Trying to add an extra element to a temporal array outside of the range."
30285 m_down->m_add_tracking(neighbor.
get_id(), p->
get_id());
30288 m->array_double_tmp[nviruses_tmp] =
30289 (1.0 - p->get_susceptibility_reduction(v, m)) *
30290 v->get_prob_infecting(m) *
30291 (1.0 - neighbor.get_transmission_reduction(v, m))
30294 m->array_virus_tmp[nviruses_tmp++] = &(*v);
30299 int which = roulette(nviruses_tmp, m);
30305 *m->array_virus_tmp[which],
30314 template<
typename TSeq>
30320 auto & v = p->get_virus();
30323 if (m->runif() < 1.0/(v->get_incubation(m)))
30328 GET_MODEL(m, model);
30339 template<
typename TSeq>
30344 GET_MODEL(m, model);
30347 bool detected = m->runif() < 1.0/m->par(
"Days undetected");
30353 model->agent_quarantine_triggered[p->
get_id()] =
30359 bool isolation_detected = (m->par(
"Isolation period") >= 0) &&
30361 (model->isolation_willingness[p->
get_id()])
30365 if (isolation_detected)
30369 const auto & v = p->get_virus();
30370 m->array_double_tmp[0] = 1.0 - (1.0 - v->get_prob_recovery(m)) *
30371 (1.0 - p->get_recovery_enhancer(v, m));
30372 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
30378 if (isolation_detected)
30393 else if (which == 1)
30410 else if ((which == 2) && isolation_detected)
30423 template<
typename TSeq>
30428 GET_MODEL(m, model);
30432 int days_since = m->
today() - model->day_onset[p->
get_id()];
30435 (m->par(
"Isolation period") <= days_since) ?
30439 m->array_double_tmp[0] = 1.0 -
30440 (1.0 - p->get_virus()->get_prob_recovery(m)) *
30441 (1.0 - p->get_recovery_enhancer(p->get_virus(), m));
30444 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
30463 else if (which == 1)
30479 else if ((which == 2) && unisolate)
30489 template<
typename TSeq>
30494 GET_MODEL(m, model);
30498 int days_since = m->
today() - model->day_flagged[p->
get_id()];
30500 bool unquarantine =
30501 (m->par(
"Quarantine period") <= days_since) ?
30513 template<
typename TSeq>
30518 GET_MODEL(m, model);
30522 int days_since = m->
today() - model->day_flagged[p->
get_id()];
30524 bool unquarantine =
30525 (m->par(
"Quarantine period") <= days_since) ?
30528 if (m->runif() < 1.0/(p->get_virus()->get_incubation(m)))
30549 else if (unquarantine)
30558 template<
typename TSeq>
30563 GET_MODEL(m, model);
30567 int days_since = m->
today() - model->day_onset[p->
get_id()];
30570 (m->par(
"Isolation period") <= days_since) ?
30582 template<
typename TSeq>
30588 if (m->runif() < 1.0/m->par(
"Hospitalization period"))
30593 template<
typename TSeq>
30597 for (
size_t agent_i = 0u; agent_i < Model<TSeq>::size(); ++agent_i)
30603 agent_quarantine_triggered[agent_i] !=
30608 if (this->par(
"Quarantine period") < 0)
30614 size_t n_contacts = this->tracking_matrix_size[agent_i];
30615 if (n_contacts >= EPI_MAX_TRACKING)
30616 n_contacts = EPI_MAX_TRACKING;
30618 auto success_rate = this->par(
"Contact tracing success rate");
30619 for (
size_t contact_i = 0u; contact_i < n_contacts; ++contact_i)
30626 size_t contact_id = this->tracking_matrix[
30632 if (agent.get_state() > INFECTED)
30636 if (agent.get_n_tools() != 0u)
30639 if (quarantine_willingness[contact_id])
30642 switch (agent.get_state())
30645 agent.change_state(
this, QUARANTINED_SUSCEPTIBLE);
30649 agent.change_state(
this, QUARANTINED_EXPOSED);
30653 if (isolation_willingness[contact_id])
30655 agent.change_state(
this, ISOLATED);
30660 throw std::logic_error(
30661 "The agent is not in a state that can be quarantined."
30669 agent_quarantine_triggered[agent_i] = QUARANTINE_PROCESS_DONE;
30697 template<
typename TSeq>
30700 const std::string & vname,
30701 epiworld_fast_uint n,
30702 epiworld_double prevalence,
30703 epiworld_double contact_rate,
30704 epiworld_double transmission_rate,
30705 epiworld_double avg_incubation_days,
30706 epiworld_double recovery_rate,
30707 std::vector< double > contact_matrix,
30708 epiworld_double hospitalization_rate,
30709 epiworld_double hospitalization_period,
30711 epiworld_double days_undetected,
30712 epiworld_fast_int quarantine_period,
30713 epiworld_double quarantine_willingness,
30714 epiworld_double isolation_willingness,
30715 epiworld_fast_int isolation_period,
30716 epiworld_double contact_tracing_success_rate,
30717 epiworld_fast_uint contact_tracing_days_prior
30722 this->contact_matrix = contact_matrix;
30725 model.add_param(contact_rate,
"Contact rate");
30726 model.add_param(transmission_rate,
"Prob. Transmission");
30727 model.add_param(recovery_rate,
"Prob. Recovery");
30728 model.add_param(avg_incubation_days,
"Avg. Incubation days");
30729 model.add_param(hospitalization_rate,
"Hospitalization rate");
30730 model.add_param(hospitalization_period,
"Hospitalization period");
30731 model.add_param(days_undetected,
"Days undetected");
30732 model.add_param(quarantine_period,
"Quarantine period");
30734 quarantine_willingness,
"Quarantine willingness"
30737 isolation_willingness,
"Isolation willingness"
30739 model.add_param(isolation_period,
"Isolation period");
30741 contact_tracing_success_rate,
"Contact tracing success rate"
30744 contact_tracing_days_prior,
"Contact tracing days prior"
30748 model.add_state(
"Susceptible", m_update_susceptible);
30749 model.add_state(
"Exposed", m_update_exposed);
30750 model.add_state(
"Infected", m_update_infected);
30751 model.add_state(
"Isolated", m_update_isolated);
30752 model.add_state(
"Detected Hospitalized", m_update_hospitalized);
30753 model.add_state(
"Quarantined Susceptible", m_update_quarantine_suscep);
30754 model.add_state(
"Quarantined Exposed", m_update_quarantine_exposed);
30755 model.add_state(
"Isolated Recovered", m_update_isolated_recovered);
30756 model.add_state(
"Hospitalized", m_update_hospitalized);
30757 model.add_state(
"Recovered");
30760 model.add_globalevent(this->m_update_model,
"Update infected individuals");
30761 model.queuing_off();
30771 virus.set_prob_infecting(&model(
"Prob. Transmission"));
30772 virus.set_prob_recovery(&model(
"Prob. Recovery"));
30773 virus.set_incubation(&model(
"Avg. Incubation days"));
30775 model.add_virus(virus);
30777 model.queuing_off();
30780 model.agents_empty_graph(n);
30782 model.set_name(
"SEIR with Mixing and Quarantine");
30788 template<
typename TSeq>
30790 const std::string & vname,
30791 epiworld_fast_uint n,
30792 epiworld_double prevalence,
30793 epiworld_double contact_rate,
30794 epiworld_double transmission_rate,
30795 epiworld_double avg_incubation_days,
30796 epiworld_double recovery_rate,
30797 std::vector< double > contact_matrix,
30798 epiworld_double hospitalization_rate,
30799 epiworld_double hospitalization_period,
30801 epiworld_double days_undetected,
30802 epiworld_fast_int quarantine_period,
30803 epiworld_double quarantine_willingness,
30804 epiworld_double isolation_willingness,
30805 epiworld_fast_int isolation_period,
30806 epiworld_double contact_tracing_success_rate,
30807 epiworld_fast_uint contact_tracing_days_prior
30811 this->contact_matrix = contact_matrix;
30820 avg_incubation_days,
30823 hospitalization_rate,
30824 hospitalization_period,
30828 quarantine_willingness,
30829 isolation_willingness,
30831 contact_tracing_success_rate,
30832 contact_tracing_days_prior
30839 template<
typename TSeq>
30841 std::vector< double > proportions_,
30847 create_init_function_seir<TSeq>(proportions_)
30855 #undef SAMPLE_FROM_PROBS
30889 [[assume((output) !=
nullptr)]];
30892 #define GET_MODEL(model, output) \
30893 auto * output = dynamic_cast< ModelMeaslesMixing<TSeq> * >( (model) ); \
30894 assert((output) != nullptr);
30897 #define SAMPLE_FROM_PROBS(n, ans) \
30899 epiworld_double p_total = m->runif(); \
30900 for (ans = 0u; ans < n; ++ans) \
30902 if (p_total < m->array_double_tmp[ans]) \
30904 m->array_double_tmp[ans + 1] += m->array_double_tmp[ans]; \
30944 template<
typename TSeq = EPI_DEFAULT_TSEQ>
30950 std::vector< size_t > infectious;
30953 std::vector< size_t > n_infectious_per_group;
30956 std::vector< size_t > entity_indices;
30958 void m_update_infectious_list();
30959 std::vector< size_t > sampled_agents;
30960 size_t sample_agents(
30962 std::vector< size_t > & sampled_agents
30964 std::vector< double > adjusted_contact_rate;
30965 std::vector< double > contact_matrix;
30968 std::vector< int > sampled_sizes;
30985 std::vector< bool > quarantine_willingness;
30986 std::vector< bool > isolation_willingness;
30987 std::vector< size_t > agent_quarantine_triggered;
30988 std::vector< int > day_flagged;
30989 std::vector< int > day_rash_onset;
30990 std::vector< int > day_exposed;
30992 void m_quarantine_process();
30996 std::vector< size_t > tracking_matrix;
30997 std::vector< size_t > tracking_matrix_size;
30998 std::vector< size_t > tracking_matrix_date;
31000 void m_add_tracking(
size_t infectious_id,
size_t agent_id);
31004 static const int SUSCEPTIBLE = 0;
31005 static const int EXPOSED = 1;
31006 static const int PRODROMAL = 2;
31007 static const int RASH = 3;
31008 static const int ISOLATED = 4;
31009 static const int ISOLATED_RECOVERED = 5;
31010 static const int DETECTED_HOSPITALIZED = 6;
31011 static const int QUARANTINED_EXPOSED = 7;
31012 static const int QUARANTINED_SUSCEPTIBLE = 8;
31013 static const int QUARANTINED_PRODROMAL = 9;
31014 static const int QUARANTINED_RECOVERED = 10;
31015 static const int HOSPITALIZED = 11;
31016 static const int RECOVERED = 12;
31018 static const size_t QUARANTINE_PROCESS_INACTIVE = 0u;
31019 static const size_t QUARANTINE_PROCESS_ACTIVE = 1u;
31020 static const size_t QUARANTINE_PROCESS_DONE = 2u;
31052 epiworld_fast_uint n,
31053 epiworld_double prevalence,
31054 epiworld_double contact_rate,
31055 epiworld_double transmission_rate,
31056 epiworld_double vax_efficacy,
31057 epiworld_double vax_reduction_recovery_rate,
31058 epiworld_double incubation_period,
31059 epiworld_double prodromal_period,
31060 epiworld_double rash_period,
31061 std::vector< double > contact_matrix,
31062 epiworld_double hospitalization_rate,
31063 epiworld_double hospitalization_period,
31065 epiworld_double days_undetected,
31066 epiworld_fast_int quarantine_period,
31067 epiworld_double quarantine_willingness,
31068 epiworld_double isolation_willingness,
31069 epiworld_fast_int isolation_period,
31070 epiworld_double prop_vaccinated,
31071 epiworld_double contact_tracing_success_rate = 1.0,
31072 epiworld_fast_uint contact_tracing_days_prior = 4u
31100 epiworld_fast_uint n,
31101 epiworld_double prevalence,
31102 epiworld_double contact_rate,
31103 epiworld_double transmission_rate,
31104 epiworld_double vax_efficacy,
31105 epiworld_double vax_reduction_recovery_rate,
31106 epiworld_double incubation_period,
31107 epiworld_double prodromal_period,
31108 epiworld_double rash_period,
31109 std::vector< double > contact_matrix,
31110 epiworld_double hospitalization_rate,
31111 epiworld_double hospitalization_period,
31113 epiworld_double days_undetected,
31114 epiworld_fast_int quarantine_period,
31115 epiworld_double quarantine_willingness,
31116 epiworld_double isolation_willingness,
31117 epiworld_fast_int isolation_period,
31118 epiworld_double prop_vaccinated,
31119 epiworld_double contact_tracing_success_rate = 1.0,
31120 epiworld_fast_uint contact_tracing_days_prior = 4u
31130 epiworld_fast_uint ndays,
31153 std::vector< double > proportions_,
31154 std::vector< int > queue_ = {}
31163 contact_matrix = cmat;
31173 return contact_matrix;
31182 return agent_quarantine_triggered;
31191 return quarantine_willingness;
31200 return isolation_willingness;
31205 template<
typename TSeq>
31208 GET_MODEL(m, model);
31209 model->m_quarantine_process();
31210 model->events_run();
31211 model->m_update_infectious_list();
31215 template<
typename TSeq>
31217 size_t infectious_id,
31224 agent_quarantine_triggered[infectious_id] >=
31230 size_t loc = tracking_matrix_size[infectious_id] % EPI_MAX_TRACKING;
31232 tracking_matrix[loc] = agent_id;
31236 tracking_matrix_size[infectious_id]++;
31242 template<
typename TSeq>
31248 std::fill(n_infectious_per_group.begin(), n_infectious_per_group.end(), 0u);
31250 for (
const auto & a : agents)
31255 if (a.get_n_entities() > 0u)
31257 const auto & entity = a.get_entity(0u);
31260 entity_indices[entity.get_id()] +
31262 n_infectious_per_group[entity.get_id()]++
31274 template<
typename TSeq>
31277 std::vector< size_t > & sampled_agents
31281 size_t agent_group_id = agent->get_entity(0u).get_id();
31282 size_t ngroups = this->entities.size();
31285 for (
size_t g = 0; g < ngroups; ++g)
31288 size_t group_size = n_infectious_per_group[g];
31290 if (group_size == 0u)
31296 adjusted_contact_rate[g] * contact_matrix[
31297 MM(agent_group_id, g, ngroups)
31305 for (
int s = 0; s < nsamples; ++s)
31312 if (which >=
static_cast<int>(group_size))
31313 which =
static_cast<int>(group_size) - 1;
31316 auto & a = this->population.at(infectious.at(entity_indices[g] + which));
31318 auto & a = this->get_agent(infectious[entity_indices[g] + which]);
31323 throw std::logic_error(
31324 "The agent is not in prodromal state, but it should be."
31329 if (a.get_id() == agent->
get_id())
31332 sampled_agents[samp_id++] = a.get_id();
31342 template<
typename TSeq>
31344 epiworld_fast_uint ndays,
31355 template<
typename TSeq>
31362 size_t nentities = this->entities.size();
31363 if (this->contact_matrix.size() != nentities*nentities)
31364 throw std::length_error(
31365 std::string(
"The contact matrix must be a square matrix of size ") +
31366 std::string(
"nentities x nentities. ") +
31367 std::to_string(this->contact_matrix.size()) +
31368 std::string(
" != ") + std::to_string(nentities*nentities) +
31372 for (
size_t i = 0u; i < this->entities.size(); ++i)
31375 for (
size_t j = 0u; j < this->entities.size(); ++j)
31377 if (this->contact_matrix[MM(i, j, nentities)] < 0.0)
31378 throw std::range_error(
31379 std::string(
"The contact matrix must be non-negative. ") +
31380 std::to_string(this->contact_matrix[MM(i, j, nentities)]) +
31381 std::string(
" < 0.")
31383 sum += this->contact_matrix[MM(i, j, nentities)];
31385 if (sum < 0.999 || sum > 1.001)
31386 throw std::range_error(
31387 std::string(
"The contact matrix must have rows that add to one. ") +
31388 std::to_string(sum) +
31389 std::string(
" != 1.")
31397 n_infectious_per_group.resize(this->entities.size(), 0u);
31398 std::fill(n_infectious_per_group.begin(), n_infectious_per_group.end(), 0u);
31402 std::fill(infectious.begin(), infectious.end(), 0u);
31405 entity_indices.resize(this->entities.size(), 0u);
31406 std::fill(entity_indices.begin(), entity_indices.end(), 0u);
31407 for (
size_t i = 1u; i < this->entities.size(); ++i)
31410 entity_indices[i] +=
31411 this->entities[i - 1].size() +
31412 entity_indices[i - 1]
31418 adjusted_contact_rate.clear();
31419 adjusted_contact_rate.resize(this->entities.size(), 0.0);
31421 for (
size_t i = 0u; i < this->entities.size(); ++i)
31424 adjusted_contact_rate[i] =
31426 static_cast< epiworld_double
> (this->get_entity(i).size());
31430 if (adjusted_contact_rate[i] > 1.0)
31431 adjusted_contact_rate[i] = 1.0;
31435 this->m_update_infectious_list();
31438 quarantine_willingness.resize(this->size(),
false);
31439 isolation_willingness.resize(this->size(),
false);
31440 for (
size_t idx = 0; idx < quarantine_willingness.size(); ++idx)
31442 quarantine_willingness[idx] =
31444 isolation_willingness[idx] =
31448 agent_quarantine_triggered.resize(this->size(), 0u);
31450 agent_quarantine_triggered.begin(),
31451 agent_quarantine_triggered.end(),
31455 day_flagged.resize(this->size(), 0);
31457 day_flagged.begin(),
31462 day_rash_onset.resize(this->size(), 0);
31464 day_rash_onset.begin(),
31465 day_rash_onset.end(),
31469 day_exposed.resize(this->size(), 0);
31471 day_exposed.begin(),
31478 std::fill(tracking_matrix.begin(), tracking_matrix.end(), 0u);
31481 std::fill(tracking_matrix_size.begin(), tracking_matrix_size.end(), 0u);
31484 std::fill(tracking_matrix_date.begin(), tracking_matrix_date.end(), 0u);
31490 template<
typename TSeq>
31498 #if __cplusplus >= 202302L
31500 [[assume(ptr !=
nullptr)]];
31503 assert(ptr !=
nullptr);
31510 template<
typename TSeq>
31515 if (p->get_n_entities() == 0)
31520 GET_MODEL(m, m_down);
31522 size_t ndraws = m_down->sample_agents(p, m_down->sampled_agents);
31525 m_down->sampled_sizes.push_back(
static_cast<int>(ndraws));
31532 int nviruses_tmp = 0;
31533 for (
size_t n = 0u; n < ndraws; ++n)
31536 auto & neighbor = m->get_agent(m_down->sampled_agents[n]);
31538 auto & v = neighbor.get_virus();
31541 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
31542 throw std::logic_error(
31543 "Trying to add an extra element to a temporal array outside of the range."
31548 m_down->m_add_tracking(neighbor.
get_id(), p->
get_id());
31551 m->array_double_tmp[nviruses_tmp] =
31552 (1.0 - p->get_susceptibility_reduction(v, m)) *
31553 v->get_prob_infecting(m) *
31554 (1.0 - neighbor.get_transmission_reduction(v, m))
31557 m->array_virus_tmp[nviruses_tmp++] = &(*v);
31562 int which = roulette(nviruses_tmp, m);
31568 *m->array_virus_tmp[which],
31577 template<
typename TSeq>
31583 auto & v = p->get_virus();
31586 if (m->runif() < 1.0/(v->get_incubation(m)))
31599 template<
typename TSeq>
31604 GET_MODEL(m, model);
31607 if (m->runif() < 1.0/m->par(
"Prodromal period"))
31617 template<
typename TSeq>
31622 GET_MODEL(m, model);
31625 bool detected =
false;
31627 (m->par(
"Isolation period") >= 0) &&
31628 (m->runif() < 1.0/m->par(
"Days undetected"))
31631 model->agent_quarantine_triggered[p->
get_id()] =
31637 m->array_double_tmp[0] = 1.0/m->par(
"Rash period");
31638 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
31651 else if (which == 1)
31660 else if (which != 0)
31662 throw std::logic_error(
"The roulette returned an unexpected value.");
31664 else if ((which == 0u) && detected)
31676 template<
typename TSeq>
31681 GET_MODEL(m, model);
31685 int days_since = m->
today() - model->day_rash_onset[p->
get_id()];
31688 (m->par(
"Isolation period") <= days_since) ?
31692 m->array_double_tmp[0] = 1.0/m->par(
"Rash period");
31695 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
31714 else if (which == 1)
31730 else if ((which == 0) && unisolate)
31740 template<
typename TSeq>
31745 GET_MODEL(m, model);
31749 int days_since = m->
today() - model->day_flagged[p->
get_id()];
31751 bool unquarantine =
31752 (m->par(
"Quarantine period") <= days_since) ?
31764 template<
typename TSeq>
31769 GET_MODEL(m, model);
31773 int days_since = m->
today() - model->day_flagged[p->
get_id()];
31775 bool unquarantine =
31776 (m->par(
"Quarantine period") <= days_since) ?
31779 if (m->runif() < 1.0/(p->get_virus()->get_incubation(m)))
31797 else if (unquarantine)
31806 template<
typename TSeq>
31811 GET_MODEL(m, model);
31815 int days_since = m->
today() - model->day_flagged[p->
get_id()];
31817 bool unquarantine =
31818 (m->par(
"Quarantine period") <= days_since) ?
31822 if (m->runif() < (1.0/m->par(
"Prodromal period")))
31837 template<
typename TSeq>
31842 GET_MODEL(m, model);
31843 int days_since = m->
today() - model->day_flagged[p->
get_id()];
31845 if (days_since >= m->par(
"Quarantine period"))
31850 template<
typename TSeq>
31855 GET_MODEL(m, model);
31859 int days_since = m->
today() - model->day_rash_onset[p->
get_id()];
31862 (m->par(
"Isolation period") <= days_since) ?
31874 template<
typename TSeq>
31880 if (m->runif() < 1.0/m->par(
"Hospitalization period"))
31885 template<
typename TSeq>
31889 for (
size_t agent_i = 0u; agent_i < Model<TSeq>::size(); ++agent_i)
31895 agent_quarantine_triggered[agent_i] !=
31900 if (this->par(
"Quarantine period") < 0)
31906 size_t n_contacts = this->tracking_matrix_size[agent_i];
31907 if (n_contacts >= EPI_MAX_TRACKING)
31908 n_contacts = EPI_MAX_TRACKING;
31911 size_t day_rash_onset_agent_i = this->day_rash_onset[agent_i];
31913 for (
size_t contact_i = 0u; contact_i < n_contacts; ++contact_i)
31918 bool within_days_prior =
31919 (day_rash_onset_agent_i - tracking_matrix_date[loc]) <=
31921 if (!within_days_prior)
31928 size_t contact_id = this->tracking_matrix[loc];
31932 if (agent.get_state() > RASH)
31936 if (agent.get_n_tools() != 0u)
31940 quarantine_willingness[contact_id] &&
31944 switch (agent.get_state())
31947 agent.change_state(
this, QUARANTINED_SUSCEPTIBLE);
31951 agent.change_state(
this, QUARANTINED_EXPOSED);
31955 agent.change_state(
this, QUARANTINED_PRODROMAL);
31959 if (isolation_willingness[contact_id])
31961 agent.change_state(
this, ISOLATED);
31966 throw std::logic_error(
31967 "The agent is not in a state that can be quarantined."
31975 agent_quarantine_triggered[agent_i] =
32007 template<
typename TSeq>
32010 epiworld_fast_uint n,
32011 epiworld_double prevalence,
32012 epiworld_double contact_rate,
32013 epiworld_double transmission_rate,
32014 epiworld_double vax_efficacy,
32015 epiworld_double vax_reduction_recovery_rate,
32016 epiworld_double incubation_period,
32017 epiworld_double prodromal_period,
32018 epiworld_double rash_period,
32019 std::vector< double > contact_matrix,
32020 epiworld_double hospitalization_rate,
32021 epiworld_double hospitalization_period,
32023 epiworld_double days_undetected,
32024 epiworld_fast_int quarantine_period,
32025 epiworld_double quarantine_willingness,
32026 epiworld_double isolation_willingness,
32027 epiworld_fast_int isolation_period,
32028 epiworld_double prop_vaccinated,
32029 epiworld_double contact_tracing_success_rate,
32030 epiworld_fast_uint contact_tracing_days_prior
32035 this->contact_matrix = contact_matrix;
32038 model.add_param(contact_rate,
"Contact rate");
32039 model.add_param(transmission_rate,
"Transmission rate");
32040 model.add_param(incubation_period,
"Incubation period");
32041 model.add_param(prodromal_period,
"Prodromal period");
32042 model.add_param(rash_period,
"Rash period");
32043 model.add_param(hospitalization_rate,
"Hospitalization rate");
32044 model.add_param(hospitalization_period,
"Hospitalization period");
32045 model.add_param(days_undetected,
"Days undetected");
32046 model.add_param(quarantine_period,
"Quarantine period");
32048 quarantine_willingness,
"Quarantine willingness"
32051 isolation_willingness,
"Isolation willingness"
32053 model.add_param(isolation_period,
"Isolation period");
32055 contact_tracing_success_rate,
"Contact tracing success rate"
32058 contact_tracing_days_prior,
"Contact tracing days prior"
32060 model.add_param(prop_vaccinated,
"Vaccination rate");
32061 model.add_param(vax_efficacy,
"Vax efficacy");
32062 model.add_param(vax_reduction_recovery_rate,
"(IGNORED) Vax improved recovery");
32065 model.add_state(
"Susceptible", m_update_susceptible);
32066 model.add_state(
"Exposed", m_update_exposed);
32067 model.add_state(
"Prodromal", m_update_prodromal);
32068 model.add_state(
"Rash", m_update_rash);
32069 model.add_state(
"Isolated", m_update_isolated);
32070 model.add_state(
"Isolated Recovered", m_update_isolated_recovered);
32071 model.add_state(
"Detected Hospitalized", m_update_hospitalized);
32072 model.add_state(
"Quarantined Exposed", m_update_quarantine_exposed);
32073 model.add_state(
"Quarantined Susceptible", m_update_quarantine_suscep);
32074 model.add_state(
"Quarantined Prodromal", m_update_quarantine_prodromal);
32075 model.add_state(
"Quarantined Recovered", m_update_quarantine_recovered);
32076 model.add_state(
"Hospitalized", m_update_hospitalized);
32077 model.add_state(
"Recovered");
32080 model.add_globalevent(this->m_update_model,
"Update infected individuals");
32081 model.queuing_off();
32091 virus.set_prob_infecting(&model(
"Transmission rate"));
32092 virus.set_prob_recovery(&model(
"Rash period"));
32093 virus.set_incubation(&model(
"Incubation period"));
32095 model.add_virus(virus);
32098 Tool<> vaccine(
"Vaccine");
32099 vaccine.set_susceptibility_reduction(&model(
"Vax efficacy"));
32100 vaccine.set_recovery_enhancer(&model(
"(IGNORED) Vax improved recovery"));
32101 vaccine.set_distribution(
32102 distribute_tool_randomly(prop_vaccinated,
true)
32105 model.add_tool(vaccine);
32107 model.queuing_off();
32110 model.agents_empty_graph(n);
32112 model.set_name(
"Measles with Mixing and Quarantine");
32118 template<
typename TSeq>
32120 epiworld_fast_uint n,
32121 epiworld_double prevalence,
32122 epiworld_double contact_rate,
32123 epiworld_double transmission_rate,
32124 epiworld_double vax_efficacy,
32125 epiworld_double vax_reduction_recovery_rate,
32126 epiworld_double incubation_period,
32127 epiworld_double prodromal_period,
32128 epiworld_double rash_period,
32129 std::vector< double > contact_matrix,
32130 epiworld_double hospitalization_rate,
32131 epiworld_double hospitalization_period,
32133 epiworld_double days_undetected,
32134 epiworld_fast_int quarantine_period,
32135 epiworld_double quarantine_willingness,
32136 epiworld_double isolation_willingness,
32137 epiworld_fast_int isolation_period,
32138 epiworld_double prop_vaccinated,
32139 epiworld_double contact_tracing_success_rate,
32140 epiworld_fast_uint contact_tracing_days_prior
32144 this->contact_matrix = contact_matrix;
32153 vax_reduction_recovery_rate,
32158 hospitalization_rate,
32159 hospitalization_period,
32163 quarantine_willingness,
32164 isolation_willingness,
32167 contact_tracing_success_rate,
32168 contact_tracing_days_prior
32175 template<
typename TSeq>
32177 std::vector< double > proportions_,
32183 create_init_function_seir<TSeq>(proportions_)
32191 #undef SAMPLE_FROM_PROBS
32223 [[assume((output) !=
nullptr)]];
32226 #define GET_MODEL(model, output) \
32227 auto * output = dynamic_cast< ModelMeaslesMixingRiskQuarantine<TSeq> * >( (model) ); \
32228 assert((output) != nullptr);
32236 #define SAMPLE_FROM_PROBS(n, ans) \
32238 epiworld_double p_total = m->runif(); \
32239 for (ans = 0; ans < static_cast<int>(n); ans++) \
32241 if (p_total < m->array_double_tmp[ans]) \
32243 m->array_double_tmp[ans + 1] += m->array_double_tmp[ans]; \
32270 template<
typename TSeq = EPI_DEFAULT_TSEQ>
32275 std::vector< size_t > infectious;
32278 std::vector< size_t > n_infectious_per_group;
32281 std::vector< size_t > infectious_entity_indices;
32283 double m_get_risk_period(
size_t agent_id);
32286 void m_update_infectious_list();
32289 std::vector< size_t > sampled_agents;
32297 size_t sample_infectious_agents(
32299 std::vector< size_t > & sampled_agents
32302 std::vector< double > adjusted_contact_rate;
32305 std::vector< double > contact_matrix;
32308 std::vector< int > sampled_sizes;
32325 std::vector< bool > quarantine_willingness;
32326 std::vector< bool > isolation_willingness;
32327 std::vector< int > day_flagged;
32328 std::vector< int > day_rash_onset;
32330 std::vector< size_t > agents_triggered_contact_tracing;
32331 size_t agents_triggered_contact_tracing_size = 0u;
32332 std::vector< bool > agents_checked_contact_tracing;
32335 void m_add_contact_tracing(
size_t agent_id);
32337 std::vector< int > quarantine_risk_level;
32339 void m_quarantine_process();
32346 std::vector< size_t > days_quarantine_triggered;
32350 static constexpr
int SUSCEPTIBLE = 0;
32351 static constexpr
int EXPOSED = 1;
32352 static constexpr
int PRODROMAL = 2;
32353 static constexpr
int RASH = 3;
32354 static constexpr
int ISOLATED = 4;
32355 static constexpr
int ISOLATED_RECOVERED = 5;
32356 static constexpr
int DETECTED_HOSPITALIZED = 6;
32357 static constexpr
int QUARANTINED_EXPOSED = 7;
32358 static constexpr
int QUARANTINED_SUSCEPTIBLE = 8;
32359 static constexpr
int QUARANTINED_PRODROMAL = 9;
32360 static constexpr
int QUARANTINED_RECOVERED = 10;
32361 static constexpr
int HOSPITALIZED = 11;
32362 static constexpr
int RECOVERED = 12;
32364 static constexpr
size_t QUARANTINE_PROCESS_INACTIVE = 0u;
32365 static constexpr
size_t QUARANTINE_PROCESS_ACTIVE = 1u;
32366 static constexpr
size_t QUARANTINE_PROCESS_DONE = 2u;
32369 static constexpr
int RISK_LOW = 0;
32370 static constexpr
int RISK_MEDIUM = 1;
32371 static constexpr
int RISK_HIGH = 2;
32405 epiworld_fast_uint n,
32406 epiworld_double prevalence,
32407 epiworld_double contact_rate,
32408 epiworld_double transmission_rate,
32409 epiworld_double vax_efficacy,
32410 epiworld_double incubation_period,
32411 epiworld_double prodromal_period,
32412 epiworld_double rash_period,
32413 std::vector< double > contact_matrix,
32414 epiworld_double hospitalization_rate,
32415 epiworld_double hospitalization_period,
32417 epiworld_double days_undetected,
32418 epiworld_fast_int quarantine_period_high,
32419 epiworld_fast_int quarantine_period_medium,
32420 epiworld_fast_int quarantine_period_low,
32421 epiworld_double quarantine_willingness,
32422 epiworld_double isolation_willingness,
32423 epiworld_fast_int isolation_period,
32424 epiworld_double prop_vaccinated,
32425 epiworld_double detection_rate_quarantine,
32426 epiworld_double contact_tracing_success_rate = 1.0,
32427 epiworld_fast_uint contact_tracing_days_prior = 4u
32457 epiworld_fast_uint n,
32458 epiworld_double prevalence,
32459 epiworld_double contact_rate,
32460 epiworld_double transmission_rate,
32461 epiworld_double vax_efficacy,
32462 epiworld_double incubation_period,
32463 epiworld_double prodromal_period,
32464 epiworld_double rash_period,
32465 std::vector< double > contact_matrix,
32466 epiworld_double hospitalization_rate,
32467 epiworld_double hospitalization_period,
32469 epiworld_double days_undetected,
32470 epiworld_fast_int quarantine_period_high,
32471 epiworld_fast_int quarantine_period_medium,
32472 epiworld_fast_int quarantine_period_low,
32473 epiworld_double quarantine_willingness,
32474 epiworld_double isolation_willingness,
32475 epiworld_fast_int isolation_period,
32476 epiworld_double prop_vaccinated,
32477 epiworld_double detection_rate_quarantine,
32478 epiworld_double contact_tracing_success_rate = 1.0,
32479 epiworld_fast_uint contact_tracing_days_prior = 4u
32489 epiworld_fast_uint ndays,
32512 std::vector< double > proportions_,
32513 std::vector< int > queue_ = {}
32522 contact_matrix = cmat;
32532 return contact_matrix;
32541 return quarantine_willingness;
32550 return isolation_willingness;
32559 return quarantine_risk_level;
32569 return days_quarantine_triggered;
32576 template<
typename TSeq>
32581 GET_MODEL(m, model);
32584 model->m_quarantine_process();
32588 model->m_update_infectious_list();
32593 template<
typename TSeq>
32597 double risk_period = 0.0;
32600 int risk_level = quarantine_risk_level[agent_id];
32602 if (risk_level == RISK_HIGH)
32603 risk_period = this->par(
"Quarantine period high");
32604 else if (risk_level == RISK_MEDIUM)
32605 risk_period = this->par(
"Quarantine period medium");
32607 risk_period = this->par(
"Quarantine period low");
32609 return risk_period;
32613 template<
typename TSeq>
32621 std::fill(n_infectious_per_group.begin(), n_infectious_per_group.end(), 0u);
32624 for (
const auto & a : agents)
32627 if (a.get_state() == PRODROMAL)
32629 if (a.get_n_entities() > 0u)
32631 const auto & entity = a.get_entity(0u);
32634 infectious_entity_indices[entity.get_id()] +
32636 n_infectious_per_group[entity.get_id()]++
32648 template<
typename TSeq>
32651 std::vector< size_t > & sampled_agents
32656 if (agent->get_n_entities() == 0u)
32659 size_t agent_group_id = agent->get_entity(0u).get_id();
32660 size_t ngroups = this->entities.size();
32663 for (
size_t g = 0; g < ngroups; ++g)
32666 size_t group_size = n_infectious_per_group[g];
32668 if (group_size == 0u)
32674 adjusted_contact_rate[g] * contact_matrix[
32675 COL_MAJOR_POS(agent_group_id, g, ngroups)
32683 for (
int s = 0; s < nsamples; ++s)
32690 if (which >=
static_cast<int>(group_size))
32691 which =
static_cast<int>(group_size) - 1;
32694 auto & a = this->population.at(
32695 infectious.at(infectious_entity_indices[g] + which)
32698 const auto & a = this->get_agent(
32699 infectious[infectious_entity_indices[g] + which]
32704 if (a.get_state() != PRODROMAL)
32705 throw std::logic_error(
32706 "The agent is not in prodromal state, but it should be."
32711 if (a.get_id() == agent->
get_id())
32714 sampled_agents[samp_id++] = a.get_id();
32724 template<
typename TSeq>
32731 if (p->get_n_entities() == 0)
32736 GET_MODEL(m, m_down);
32739 size_t ndraws = m_down->sample_infectious_agents(p, m_down->sampled_agents);
32742 m_down->sampled_sizes.push_back(
static_cast<int>(ndraws));
32749 int nviruses_tmp = 0;
32750 for (
size_t n = 0u; n < ndraws; ++n)
32754 auto & neighbor = m->get_agent(m_down->sampled_agents[n]);
32755 auto & v = neighbor.get_virus();
32758 if (nviruses_tmp >=
static_cast<int>(m->array_virus_tmp.size()))
32759 throw std::logic_error(
32760 "Trying to add an extra element to a temporal array outside of the range."
32767 m_down->contact_tracing.add_contact(neighbor.
get_id(), p->
get_id(), m->
today());
32770 m->array_double_tmp[nviruses_tmp] =
32771 (1.0 - p->get_susceptibility_reduction(v, m)) *
32772 v->get_prob_infecting(m) *
32773 (1.0 - neighbor.get_transmission_reduction(v, m))
32776 m->array_virus_tmp[nviruses_tmp++] = &(*v);
32781 int which = roulette(nviruses_tmp, m);
32786 p->set_virus(*m->array_virus_tmp[which], m, EXPOSED);
32792 template<
typename TSeq>
32798 auto & v = p->get_virus();
32801 if (m->runif() < 1.0/(v->get_incubation(m)))
32803 p->change_state(m, PRODROMAL);
32808 template<
typename TSeq>
32813 GET_MODEL(m, model);
32816 if (m->runif() < 1.0/m->par(
"Prodromal period"))
32819 bool detect_it = (model->get_days_quarantine_triggered().size() > 0u) &&
32820 (m->runif() < m->par(
"Detection rate quarantine"));
32823 p->change_state(m, detect_it ? ISOLATED : RASH);
32826 model->m_add_contact_tracing(p->
get_id());
32834 template<
typename TSeq>
32839 GET_MODEL(m, model);
32842 bool detected =
false;
32844 (m->par(
"Isolation period") >= 0) &&
32845 (m->runif() < 1.0/m->par(
"Days undetected"))
32848 model->m_add_contact_tracing(p->
get_id());
32854 m->array_double_tmp[0] = 1.0/m->par(
"Rash period");
32855 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
32861 p->rm_virus(m, detected ? ISOLATED_RECOVERED: RECOVERED);
32863 else if (which == 1)
32865 p->change_state(m, detected ? DETECTED_HOSPITALIZED : HOSPITALIZED);
32867 else if ((which == 0) && detected)
32869 p->change_state(m, ISOLATED);
32871 else if (which != 0)
32873 throw std::logic_error(
"The roulette returned an unexpected value.");
32880 template<
typename TSeq>
32885 GET_MODEL(m, model);
32889 int days_since = m->
today() - model->day_rash_onset[p->
get_id()];
32892 (m->par(
"Isolation period") <= days_since) ?
32897 m->array_double_tmp[0] = 1.0/m->par(
"Rash period");
32898 m->array_double_tmp[1] = m->par(
"Hospitalization rate");
32906 p->rm_virus(m, unisolate ? RECOVERED : ISOLATED_RECOVERED);
32909 else if (which == 1u)
32911 p->change_state(m, DETECTED_HOSPITALIZED);
32914 else if ((which == 0u) && unisolate)
32916 p->change_state(m, RASH);
32922 template<
typename TSeq>
32928 GET_MODEL(m, model);
32932 int days_since = m->
today() - model->day_rash_onset[p->
get_id()];
32934 if (m->par(
"Isolation period") <= days_since)
32935 p->change_state(m, RECOVERED);
32939 template<
typename TSeq>
32945 GET_MODEL(m, model);
32948 auto quarantine_period = model->m_get_risk_period(p->
get_id());
32951 int days_since = m->
today() - model->day_flagged[p->
get_id()];
32953 if (quarantine_period <= days_since)
32954 p->change_state(m, SUSCEPTIBLE);
32958 template<
typename TSeq>
32964 GET_MODEL(m, model);
32967 auto quarantine_period = model->m_get_risk_period(p->
get_id());
32970 int days_since = m->
today() - model->day_flagged[p->
get_id()];
32973 if (m->runif() < 1.0/(p->get_virus()->get_incubation(m)))
32979 (quarantine_period <= days_since) ?
32980 PRODROMAL : QUARANTINED_PRODROMAL
32985 else if (quarantine_period <= days_since)
32987 p->change_state(m, EXPOSED);
32992 template<
typename TSeq>
32998 GET_MODEL(m, model);
33001 auto quarantine_period = model->m_get_risk_period(p->
get_id());
33004 int days_since = m->
today() - model->day_flagged[p->
get_id()];
33007 if (m->runif() < (1.0/m->par(
"Prodromal period")))
33011 model->m_add_contact_tracing(p->
get_id());
33014 p->change_state(m, ISOLATED);
33017 else if (quarantine_period <= days_since)
33019 p->change_state(m, PRODROMAL);
33025 template<
typename TSeq>
33031 GET_MODEL(m, model);
33034 auto quarantine_period = model->m_get_risk_period(p->
get_id());
33037 int days_since = m->
today() - model->day_flagged[p->
get_id()];
33039 if (quarantine_period <= days_since)
33040 p->change_state(m, RECOVERED);
33044 template<
typename TSeq>
33051 if (m->runif() < 1.0/m->par(
"Hospitalization period"))
33052 p->rm_virus(m, RECOVERED);
33056 template<
typename TSeq>
33059 this->agents_triggered_contact_tracing[
33060 this->agents_triggered_contact_tracing_size++
33067 template<
typename TSeq>
33071 if (agents_triggered_contact_tracing_size == 0u)
33079 std::vector< bool > can_quarantine(
33088 auto state = agent.get_state();
33091 (state == SUSCEPTIBLE) ||
33092 (state == EXPOSED) ||
33093 (state == PRODROMAL) ||
33097 quarantine_risk_level[agent.
get_id()] = RISK_LOW;
33098 can_quarantine[agent.
get_id()] =
true;
33104 std::map< std::pair<size_t, size_t>,
bool> success_contact;
33109 for (
size_t i = 0u; i < agents_triggered_contact_tracing_size; ++i)
33112 size_t agent_i_idx = agents_triggered_contact_tracing[i];
33118 if (agent_i.get_n_entities() != 0u)
33120 for (
auto & agent_j_idx: agent_i.get_entity(0))
33126 agent_j.get_entity(0u).get_id() !=
33127 agent_i.get_entity(0u).get_id()
33129 throw std::logic_error(
33130 "An agent in a group has a different group id."
33136 if (can_quarantine[agent_j_idx])
33137 quarantine_risk_level[agent_j_idx] = RISK_HIGH;
33144 size_t n_contacts = contact_tracing.get_n_contacts(agent_i_idx);
33145 if (n_contacts > EPI_MAX_TRACKING)
33146 n_contacts = EPI_MAX_TRACKING;
33148 for (
size_t contact_j_idx = 0u; contact_j_idx < n_contacts; ++contact_j_idx)
33152 auto [contact_id, contact_date] = contact_tracing.get_contact(
33153 agent_i_idx, contact_j_idx
33157 if (!can_quarantine[contact_id])
33169 double days_since_contact =
33170 static_cast<double>(this->day_rash_onset[agent_i_idx]) -
33171 static_cast<double>(contact_date);
33174 if (days_since_contact > param_days_prior)
33178 auto pair = std::make_pair(agent_i_idx, contact_id);
33179 success_contact[pair] =
true;
33183 if (quarantine_risk_level[contact_id] < RISK_HIGH)
33184 quarantine_risk_level[contact_id] = RISK_MEDIUM;
33195 if (!can_quarantine[agent.
get_id()])
33199 if (agent.get_n_tools() != 0u)
33202 if (quarantine_willingness[agent.
get_id()] ==
false)
33208 auto state = agent.get_state();
33209 if (state == SUSCEPTIBLE)
33210 agent.change_state(
this, QUARANTINED_SUSCEPTIBLE);
33211 else if (state == EXPOSED)
33212 agent.change_state(
this, QUARANTINED_EXPOSED);
33213 else if (state == PRODROMAL)
33214 agent.change_state(
this, QUARANTINED_PRODROMAL);
33215 else if (state == RASH)
33217 if (isolation_willingness[agent.
get_id()])
33218 agent.change_state(
this, ISOLATED);
33221 throw std::logic_error(
33222 "An agent in an unexpected state is being quarantined."
33229 std::set< size_t > groups_ids;
33230 std::set< size_t > contacted_agents;
33231 for (
size_t i = 0u; i < agents_triggered_contact_tracing_size; ++i)
33234 auto agent_i_idx = agents_triggered_contact_tracing[i];
33239 size_t n_contacts = contact_tracing.get_n_contacts(agent_i_idx);
33240 if (n_contacts >= EPI_MAX_TRACKING)
33241 n_contacts = EPI_MAX_TRACKING;
33243 for (
size_t contact_j_idx = 0u; contact_j_idx < n_contacts; ++contact_j_idx)
33247 auto [contact_id, contact_date] = contact_tracing.get_contact(
33248 agent_i_idx, contact_j_idx
33252 auto pair = std::make_pair(agent_i_idx, contact_id);
33253 if (!success_contact[pair])
33256 contacted_agents.insert(contact_id);
33259 if (agent_i.get_n_entities() == 0u)
33262 groups_ids.insert(agent_i.get_entity(0u).get_id());
33269 for (
auto & group: groups_ids)
33275 auto state = agent_i.get_state();
33279 (state != SUSCEPTIBLE) &&
33280 (state != EXPOSED) &&
33281 (state != PRODROMAL) &&
33287 if (agent_i.get_n_tools() != 0u)
33290 if (quarantine_risk_level[agent_i_idx] != RISK_HIGH)
33291 throw std::logic_error(
33292 "An agent in a group that triggered contact tracing is not in RISK_HIGH."
33301 for (
auto & agent_i_idx: contacted_agents)
33304 auto state = agent_i.get_state();
33307 if (agent_i.get_n_tools() != 0u)
33312 if (agent_i.get_n_entities() != 0u)
33314 size_t group_id = agent_i.get_entity(0u).get_id();
33315 if (groups_ids.find(group_id) != groups_ids.end())
33319 if (quarantine_risk_level[agent_i_idx] != RISK_MEDIUM)
33320 throw std::logic_error(
33321 "An agent who was contacted by an infectious agent is not in at least RISK_MEDIUM."
33327 for (
auto i: quarantine_risk_level)
33328 n_agents_per_group[i]++;
33335 agents_triggered_contact_tracing_size = 0u;
33340 template<
typename TSeq>
33343 epiworld_fast_uint n,
33344 epiworld_double prevalence,
33345 epiworld_double contact_rate,
33346 epiworld_double transmission_rate,
33347 epiworld_double vax_efficacy,
33348 epiworld_double incubation_period,
33349 epiworld_double prodromal_period,
33350 epiworld_double rash_period,
33351 std::vector< double > contact_matrix,
33352 epiworld_double hospitalization_rate,
33353 epiworld_double hospitalization_period,
33355 epiworld_double days_undetected,
33356 epiworld_fast_int quarantine_period_high,
33357 epiworld_fast_int quarantine_period_medium,
33358 epiworld_fast_int quarantine_period_low,
33359 epiworld_double quarantine_willingness,
33360 epiworld_double isolation_willingness,
33361 epiworld_fast_int isolation_period,
33362 epiworld_double prop_vaccinated,
33363 epiworld_double detection_rate_quarantine,
33364 epiworld_double contact_tracing_success_rate,
33365 epiworld_fast_uint contact_tracing_days_prior
33370 this->contact_matrix = contact_matrix;
33373 model.add_param(contact_rate,
"Contact rate");
33374 model.add_param(transmission_rate,
"Transmission rate");
33375 model.add_param(incubation_period,
"Incubation period");
33376 model.add_param(prodromal_period,
"Prodromal period");
33377 model.add_param(rash_period,
"Rash period");
33378 model.add_param(hospitalization_rate,
"Hospitalization rate");
33379 model.add_param(hospitalization_period,
"Hospitalization period");
33380 model.add_param(days_undetected,
"Days undetected");
33381 model.add_param(quarantine_period_high,
"Quarantine period high");
33382 model.add_param(quarantine_period_medium,
"Quarantine period medium");
33383 model.add_param(quarantine_period_low,
"Quarantine period low");
33385 quarantine_willingness,
"Quarantine willingness"
33388 isolation_willingness,
"Isolation willingness"
33390 model.add_param(isolation_period,
"Isolation period");
33391 model.add_param(detection_rate_quarantine,
"Detection rate quarantine");
33393 contact_tracing_success_rate,
"Contact tracing success rate"
33396 contact_tracing_days_prior,
"Contact tracing days prior"
33398 model.add_param(prop_vaccinated,
"Vaccination rate");
33399 model.add_param(vax_efficacy,
"Vax efficacy");
33402 model.add_state(
"Susceptible", m_update_susceptible);
33403 model.add_state(
"Exposed", m_update_exposed);
33404 model.add_state(
"Prodromal", m_update_prodromal);
33405 model.add_state(
"Rash", m_update_rash);
33406 model.add_state(
"Isolated", m_update_isolated);
33407 model.add_state(
"Isolated Recovered", m_update_isolated_recovered);
33408 model.add_state(
"Detected Hospitalized", m_update_hospitalized);
33409 model.add_state(
"Quarantined Exposed", m_update_quarantine_exposed);
33410 model.add_state(
"Quarantined Susceptible", m_update_quarantine_suscep);
33411 model.add_state(
"Quarantined Prodromal", m_update_quarantine_prodromal);
33412 model.add_state(
"Quarantined Recovered", m_update_quarantine_recovered);
33413 model.add_state(
"Hospitalized", m_update_hospitalized);
33414 model.add_state(
"Recovered");
33417 model.add_globalevent(this->m_update_model,
"Update infected individuals");
33418 model.queuing_off();
33428 virus.set_prob_infecting(&model(
"Transmission rate"));
33429 virus.set_prob_recovery(&model(
"Rash period"));
33430 virus.set_incubation(&model(
"Incubation period"));
33432 model.add_virus(virus);
33435 Tool<> vaccine(
"Vaccine");
33436 vaccine.set_susceptibility_reduction(&model(
"Vax efficacy"));
33437 vaccine.set_distribution(
33438 distribute_tool_randomly(prop_vaccinated,
true)
33441 model.add_tool(vaccine);
33443 model.queuing_off();
33446 model.agents_empty_graph(n);
33448 model.set_name(
"Measles with Mixing and Risk-based Quarantine");
33454 template<
typename TSeq>
33456 epiworld_fast_uint n,
33457 epiworld_double prevalence,
33458 epiworld_double contact_rate,
33459 epiworld_double transmission_rate,
33460 epiworld_double vax_efficacy,
33461 epiworld_double incubation_period,
33462 epiworld_double prodromal_period,
33463 epiworld_double rash_period,
33464 std::vector< double > contact_matrix,
33465 epiworld_double hospitalization_rate,
33466 epiworld_double hospitalization_period,
33468 epiworld_double days_undetected,
33469 epiworld_fast_int quarantine_period_high,
33470 epiworld_fast_int quarantine_period_medium,
33471 epiworld_fast_int quarantine_period_low,
33472 epiworld_double quarantine_willingness,
33473 epiworld_double isolation_willingness,
33474 epiworld_fast_int isolation_period,
33475 epiworld_double prop_vaccinated,
33476 epiworld_double detection_rate_quarantine,
33477 epiworld_double contact_tracing_success_rate,
33478 epiworld_fast_uint contact_tracing_days_prior
33490 hospitalization_rate,
33491 hospitalization_period,
33494 quarantine_period_high,
33495 quarantine_period_medium,
33496 quarantine_period_low,
33497 quarantine_willingness,
33498 isolation_willingness,
33501 detection_rate_quarantine,
33502 contact_tracing_success_rate,
33503 contact_tracing_days_prior
33506 template<
typename TSeq>
33508 epiworld_fast_uint ndays,
33516 template<
typename TSeq>
33523 size_t nentities = this->entities.size();
33524 if (this->contact_matrix.size() != nentities*nentities)
33525 throw std::length_error(
33526 std::string(
"The contact matrix must be a square matrix of size ") +
33527 std::string(
"nentities x nentities. ") +
33528 std::to_string(this->contact_matrix.size()) +
33529 std::string(
" != ") + std::to_string(nentities*nentities) +
33533 for (
size_t i = 0u; i < this->entities.size(); ++i)
33536 for (
size_t j = 0u; j < this->entities.size(); ++j)
33538 if (this->contact_matrix[COL_MAJOR_POS(i, j, nentities)] < 0.0)
33539 throw std::range_error(
33540 std::string(
"The contact matrix must be non-negative. ") +
33541 std::to_string(this->contact_matrix[COL_MAJOR_POS(i, j, nentities)]) +
33542 std::string(
" < 0.")
33544 sum += this->contact_matrix[COL_MAJOR_POS(i, j, nentities)];
33546 if (sum < 0.999 || sum > 1.001)
33547 throw std::range_error(
33548 std::string(
"The contact matrix must have rows that add to one. ") +
33549 std::to_string(sum) +
33550 std::string(
" != 1.")
33558 n_infectious_per_group.assign(this->entities.size(), 0u);
33564 infectious_entity_indices.assign(this->entities.size(), 0u);
33565 for (
size_t i = 1u; i < this->entities.size(); ++i)
33568 infectious_entity_indices[i] +=
33569 this->entities[i - 1].size() +
33570 infectious_entity_indices[i - 1]
33576 adjusted_contact_rate.assign(this->entities.size(), 0.0);
33578 for (
size_t i = 0u; i < this->entities.size(); ++i)
33581 adjusted_contact_rate[i] =
33583 static_cast< epiworld_double
> (this->get_entity(i).size());
33587 if (adjusted_contact_rate[i] > 1.0)
33588 adjusted_contact_rate[i] = 1.0;
33592 this->m_update_infectious_list();
33595 quarantine_willingness.resize(this->size(),
false);
33596 isolation_willingness.resize(this->size(),
false);
33597 for (
size_t idx = 0; idx < quarantine_willingness.size(); ++idx)
33599 quarantine_willingness[idx] =
33601 isolation_willingness[idx] =
33605 day_flagged.assign(this->size(), 0);
33606 day_rash_onset.assign(this->size(), 0);
33610 agents_triggered_contact_tracing_size = 0u;
33611 agents_checked_contact_tracing.assign(this->size(),
false);
33613 quarantine_risk_level.assign(this->size(), RISK_LOW);
33616 contact_tracing.reset(
33622 days_quarantine_triggered.clear();
33628 template<
typename TSeq>
33640 template<
typename TSeq>
33642 std::vector< double > proportions_,
33648 create_init_function_seir<TSeq>(proportions_)
33655 #undef SAMPLE_FROM_PROBS
Definition: adjlist-bones.hpp:4
size_t vcount() const
Number of vertices/nodes in the network.
Definition: adjlist-meat.hpp:203
void read_edgelist(std::string fn, int size, int skip=0, bool directed=true)
Read an edgelist.
Definition: adjlist-meat.hpp:92
Agent (agents)
Definition: agent-bones.hpp:66
double & operator()(size_t j)
Access the j-th column of the agent.
Definition: agent-meat.hpp:836
int get_id() const
Id of the individual.
Definition: agent-meat.hpp:465
void rm_agent_by_virus(Model< TSeq > *model)
Agent removed by virus.
Definition: agent-meat.hpp:417
void swap_neighbors(Agent< TSeq > &other, size_t n_this, size_t n_other)
Swaps neighbors between the current agent and agent other
Definition: agent-meat.hpp:584
Statistical data about the process.
Definition: database-bones.hpp:34
Set of Entities (const) (useful for iterators)
Definition: entities-bones.hpp:116
Set of Entities (useful for building iterators)
Definition: entities-bones.hpp:17
Definition: entity-bones.hpp:20
Template for a Global Event.
Definition: globalevent-bones.hpp:15
Template for a Network Diffusion Model.
Definition: diffnet.hpp:15
Core class of epiworld.
Definition: model-bones.hpp:91
std::vector< Agent< TSeq > > & get_agents()
Returns a reference to the vector of agents.
Definition: model-meat.hpp:582
void run_multiple(epiworld_fast_uint ndays, epiworld_fast_uint nexperiments, int seed_=-1, std::function< void(size_t, Model< TSeq > *)> fun=make_save_run< TSeq >(), bool reset=true, bool verbose=true, int nthreads=1)
Definition: model-meat.hpp:1477
void load_agents_entities_ties(std::string fn, int skip)
Associate agents-entities from a file.
Definition: model-meat.hpp:1088
size_t get_n_tools() const
Number of tools in the model.
Definition: model-meat.hpp:1767
void write_data(std::string fn_virus_info, std::string fn_virus_hist, std::string fn_tool_info, std::string fn_tool_hist, std::string fn_total_hist, std::string fn_transmission, std::string fn_transition, std::string fn_reproductive_number, std::string fn_generation_time) const
Wrapper of DataBase::write_data
Definition: model-meat.hpp:1843
void events_run()
Executes the stored action.
Definition: model-meat.hpp:219
virtual void reset()
Reset the model.
Definition: model-meat.hpp:1962
int today() const
The current time of the model.
Definition: model-meat.hpp:1316
virtual Model< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Runs the simulation (after initialization)
Definition: model-meat.hpp:1366
size_t get_n_viruses() const
Number of viruses in the model.
Definition: model-meat.hpp:1762
Measles model with population mixing, quarantine, and contact tracing.
Definition: measlesmixing.hpp:72
ModelMeaslesMixing< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: measlesmixing.hpp:1302
ModelMeaslesMixing< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Run the model simulation.
Definition: measlesmixing.hpp:469
Model< TSeq > * clone_ptr()
Create a clone of this model.
Definition: measlesmixing.hpp:617
void reset()
Reset the model to initial state.
Definition: measlesmixing.hpp:482
Measles model with population mixing and risk-based quarantine strategies.
Definition: measlesmixingriskquarantine.hpp:62
void reset()
Reset the model to initial state.
Definition: measlesmixingriskquarantine.hpp:1307
ModelMeaslesMixingRiskQuarantine< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Run the model simulation.
Definition: measlesmixingriskquarantine.hpp:1297
ModelMeaslesMixingRiskQuarantine< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: measlesmixingriskquarantine.hpp:1431
Model< TSeq > * clone_ptr()
Create a clone of this model.
Definition: measlesmixingriskquarantine.hpp:1419
Template for a Measles model with quarantine.
Definition: measlesschool.hpp:49
void reset()
Reset the model.
Definition: seirconnected.hpp:114
Model< TSeq > * clone_ptr()
Advanced usage: Makes a copy of data and returns it as undeleted pointer.
Definition: seirconnected.hpp:125
ModelSEIRCONN< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: seirconnected.hpp:387
ModelSEIRCONN< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Runs the simulation (after initialization)
Definition: seirconnected.hpp:101
Definition: seirdconnected.hpp:6
ModelSEIRDCONN< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set up the initial states of the model.
Definition: seirdconnected.hpp:405
Definition: seirmixing.hpp:27
void reset()
Reset the model.
Definition: seirmixing.hpp:248
ModelSEIRMixing< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: seirmixing.hpp:596
ModelSEIRMixing< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Runs the simulation (after initialization)
Definition: seirmixing.hpp:235
Model< TSeq > * clone_ptr()
Advanced usage: Makes a copy of data and returns it as undeleted pointer.
Definition: seirmixing.hpp:334
SEIR model with mixing, quarantine, and contact tracing.
Definition: seirmixingquarantine.hpp:66
ModelSEIRMixingQuarantine< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: seirmixingquarantine.hpp:1205
Model< TSeq > * clone_ptr()
Create a clone of this model.
Definition: seirmixingquarantine.hpp:593
ModelSEIRMixingQuarantine< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Run the model simulation.
Definition: seirmixingquarantine.hpp:448
void reset()
Reset the model to initial state.
Definition: seirmixingquarantine.hpp:461
ModelSIRCONN< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Runs the simulation (after initialization)
Definition: sirconnected.hpp:107
Model< TSeq > * clone_ptr()
Advanced usage: Makes a copy of data and returns it as undeleted pointer.
Definition: sirconnected.hpp:131
ModelSIRCONN< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: sirconnected.hpp:360
void reset()
Reset the model.
Definition: sirconnected.hpp:119
Definition: sirdconnected.hpp:6
Definition: sirmixing.hpp:24
ModelSIRMixing< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Runs the simulation (after initialization)
Definition: sirmixing.hpp:232
ModelSIRMixing< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: sirmixing.hpp:551
Model< TSeq > * clone_ptr()
Advanced usage: Makes a copy of data and returns it as undeleted pointer.
Definition: sirmixing.hpp:324
void reset()
Reset the model.
Definition: sirmixing.hpp:244
Definition: surveillance.hpp:5
std::vector< epiworld_double > days_latent_and_infectious
Vector of days spent in latent and infectious states A row-major matrix.
Definition: surveillance.hpp:24
Controls which agents are verified at each step.
Definition: queue-bones.hpp:16
Personalized data by the user.
Definition: userdata-bones.hpp:17
Virus.
Definition: virus-bones.hpp:41
Definition: epiworld.hpp:7034
Agent (agents)
Definition: epiworld.hpp:10908
int get_id() const
Id of the individual.
Definition: epiworld.hpp:22122
void rm_agent_by_virus(Model< TSeq > *model)
Agent removed by virus.
Definition: epiworld.hpp:22074
Statistical data about the process.
Definition: epiworld.hpp:3592
Definition: epiworld.hpp:17122
Template for a Global Event.
Definition: epiworld.hpp:8148
Core class of epiworld.
Definition: epiworld.hpp:8392
std::vector< Agent< TSeq > > & get_agents()
Returns a reference to the vector of agents.
Definition: epiworld.hpp:11714
void load_agents_entities_ties(const std::vector< int > &agents_ids, const std::vector< int > &entities_ids)
Associate agents-entities from data.
Definition: epiworld.hpp:12278
void rm_globalevent(size_t i)
Remove a global action by index.
Definition: epiworld.hpp:13824
std::vector< Viruses_const< TSeq > > get_agents_viruses() const
Returns a const vector with the viruses of the agents.
Definition: epiworld.hpp:11736
void set_agents_data(double *data_, size_t ncols_)
Set the agents data object.
Definition: epiworld.hpp:13909
void queuing_on()
Activates the queuing system (default.)
Definition: epiworld.hpp:13849
bool is_queuing_on() const
Query if the queuing system is on.
Definition: epiworld.hpp:13862
void run_multiple(epiworld_fast_uint ndays, epiworld_fast_uint nexperiments, int seed_=-1, std::function< void(size_t, Model< TSeq > *)> fun=make_save_run< TSeq >(), bool reset=true, bool verbose=true, int nthreads=1)
Definition: epiworld.hpp:12609
Queue< TSeq > & get_queue()
Retrieve the Queue object.
Definition: epiworld.hpp:13868
void set_name(std::string name)
Set the name object.
Definition: epiworld.hpp:13927
void load_agents_entities_ties(std::string fn, int skip)
Associate agents-entities from a file.
Definition: epiworld.hpp:12220
size_t get_n_tools() const
Number of tools in the model.
Definition: epiworld.hpp:12899
void add_globalevent(std::function< void(Model< TSeq > *)> fun, std::string name="A global action", int date=-99)
Set a global action.
Definition: epiworld.hpp:13750
void write_data(std::string fn_virus_info, std::string fn_virus_hist, std::string fn_tool_info, std::string fn_tool_hist, std::string fn_total_hist, std::string fn_transmission, std::string fn_transition, std::string fn_reproductive_number, std::string fn_generation_time) const
Wrapper of DataBase::write_data
Definition: epiworld.hpp:12975
void events_run()
Executes the stored action.
Definition: epiworld.hpp:11351
Model< TSeq > & queuing_off()
Deactivates the queuing system.
Definition: epiworld.hpp:13855
virtual void reset()
Reset the model.
Definition: epiworld.hpp:13094
std::vector< epiworld_fast_uint > get_agents_states() const
Returns a vector with the states of the agents.
Definition: epiworld.hpp:11726
void draw(DiagramType diagram_type=DiagramType::Mermaid, const std::string &fn_output="", bool self=false)
Draws a mermaid diagram of the model.
Definition: epiworld.hpp:14108
GlobalEvent< TSeq > & get_globalevent(size_t i)
Retrieve a global action by index.
Definition: epiworld.hpp:13790
int today() const
The current time of the model.
Definition: epiworld.hpp:12448
virtual Model< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Runs the simulation (after initialization)
Definition: epiworld.hpp:12498
size_t get_n_viruses() const
Number of viruses in the model.
Definition: epiworld.hpp:12894
void rm_globalevent(std::string name)
Remove a global action by name.
Definition: epiworld.hpp:13804
std::vector< Viruses< TSeq > > get_agents_viruses()
Returns a vector with the viruses of the agents.
Definition: epiworld.hpp:11749
void set_user_data(std::vector< std::string > names)
[@
Definition: epiworld.hpp:13726
GlobalEvent< TSeq > & get_globalevent(std::string name)
Retrieve a global action by name.
Definition: epiworld.hpp:13776
void events_add(Agent< TSeq > *agent_, VirusPtr< TSeq > virus_, ToolPtr< TSeq > tool_, Entity< TSeq > *entity_, epiworld_fast_int new_state_, epiworld_fast_int queue_, EventFun< TSeq > call_, int idx_agent_, int idx_object_)
Construct a new Event object.
Definition: epiworld.hpp:11300
virtual Model< TSeq > * clone_ptr()
Advanced usage: Makes a copy of data and returns it as undeleted pointer.
Definition: epiworld.hpp:11512
A simple progress bar.
Definition: epiworld.hpp:866
Controls which agents are verified at each step.
Definition: epiworld.hpp:8002
Personalized data by the user.
Definition: epiworld.hpp:3108
Virus.
Definition: epiworld.hpp:10657
Measles model with population mixing, quarantine, and contact tracing.
Definition: epiworld.hpp:30946
void set_contact_matrix(std::vector< double > cmat)
Set the contact matrix for population mixing.
Definition: epiworld.hpp:31161
std::vector< bool > get_isolation_willingness() const
Get the isolation willingness for all agents.
Definition: epiworld.hpp:31198
std::vector< double > get_contact_matrix() const
Get the current contact matrix.
Definition: epiworld.hpp:31171
std::vector< bool > get_quarantine_willingness() const
Get the quarantine willingness for all agents.
Definition: epiworld.hpp:31189
std::vector< size_t > get_agent_quarantine_triggered() const
Get the quarantine trigger status for all agents.
Definition: epiworld.hpp:31180
Measles model with population mixing and risk-based quarantine strategies.
Definition: epiworld.hpp:32272
auto get_quarantine_willingness() const
Get the quarantine willingness for all agents.
Definition: epiworld.hpp:32539
void set_contact_matrix(std::vector< double > cmat)
Set the contact matrix for population mixing.
Definition: epiworld.hpp:32520
void reset()
Reset the model to initial state.
Definition: epiworld.hpp:33517
ModelMeaslesMixingRiskQuarantine< TSeq > & run(epiworld_fast_uint ndays, int seed=-1)
Run the model simulation.
Definition: epiworld.hpp:33507
ModelMeaslesMixingRiskQuarantine< TSeq > & initial_states(std::vector< double > proportions_, std::vector< int > queue_={})
Set the initial states of the model.
Definition: epiworld.hpp:33641
const auto & get_quarantine_risk_level() const
Get the risk level assigned to each agent for quarantine purposes.
Definition: epiworld.hpp:32557
auto get_days_quarantine_triggered() const
Get the total number of quarantines that have occurred.
Definition: epiworld.hpp:32567
Model< TSeq > * clone_ptr()
Create a clone of this model.
Definition: epiworld.hpp:33629
auto get_contact_matrix() const
Get the current contact matrix.
Definition: epiworld.hpp:32530
auto get_isolation_willingness() const
Get the isolation willingness for all agents.
Definition: epiworld.hpp:32548
Template for a Measles model with quarantine.
Definition: epiworld.hpp:28820
Definition: epiworld.hpp:25179
Definition: epiworld.hpp:26526
Definition: epiworld.hpp:27585
SEIR model with mixing, quarantine, and contact tracing.
Definition: epiworld.hpp:29701
std::vector< double > get_contact_matrix() const
Get the current contact matrix.
Definition: epiworld.hpp:29907
std::vector< bool > get_isolation_willingness() const
Get the isolation willingness for all agents.
Definition: epiworld.hpp:29934
std::vector< bool > get_quarantine_willingness() const
Get the quarantine willingness for all agents.
Definition: epiworld.hpp:29925
std::vector< size_t > get_agent_quarantine_triggered() const
Get the quarantine trigger status for all agents.
Definition: epiworld.hpp:29916
void set_contact_matrix(std::vector< double > cmat)
Set the contact matrix for population mixing.
Definition: epiworld.hpp:29897
Definition: epiworld.hpp:24738
Definition: epiworld.hpp:26169
Template for a Susceptible-Infected-Removed (SIR) model.
Definition: epiworld.hpp:26992
Definition: epiworld.hpp:28211
#define SAMPLE_FROM_PROBS(n, ans)
Macro to sample from a list of probabilities.
Definition: measlesmixingriskquarantine.hpp:26
Virus< TSeq > * sample_virus_single(Agent< TSeq > *p, Model< TSeq > *m)
Sample from neighbors pool of viruses (at most one)
Definition: epiworld.hpp:18274
std::function< Virus< TSeq > *(Agent< TSeq > *, Model< TSeq > *)> make_sample_virus_neighbors(std::vector< epiworld_fast_uint > exclude={})
Make a function to sample from neighbors.
Definition: epiworld.hpp:18106
Functions for sampling viruses.
Definition: agent-meat-virus-sampling.hpp:8